      SUBROUTINE INPUT(TITLE)  
C  
C **  SUBROUTINE INPUT READS ALL INPUT DATA EXCEPT DATA IN LXLY.INP,  
C **  MASK.INP AND RESTART.INP  
C  CHANGE RECORD  
C  ADDED BODY FORCE SWITCH ISBODYF (C14) AND READ OF INPUT FILE  
C  FBODY.INP CONTAINING THE BODY FORCE FBODYFX AND FBODYFY  
C  MODIFIED BED MECHANICS COEFFICIENT SET ON (C38)  
C  MODIFIED BEDLOAD FUNCTIONAL RELATIONSHIP TO MORE GENERAL FORM (C42A)  
C  ADDED BEDLOAD OUTFLOW/RECIRCULATION BOUNDARY CONDTION SWITCH (C42A)  
C  AND BED LOAD BOUNDARY CONDITION FILE SEDBLBC.INP  
C  ADDED SWITCH IBMECHK ON C38 TO CONTROL FORM OF FUNCTION FHYDCN  
C  REPLACED LSBDLDBC WITH LSBLBCU AND LSBLBCD AND MODIFIED FORMAT OF  
C  SEDBLBC.INP  
C  FIXED THE SPECIFICATION OF LSBLBCD FOR NO RECIRCULATION CASE  
C  INITIALIZED NSBDLDBC  
C  ADDED ICK2COR,CK2UUM,CK2VVM,CK2UVM,CK2UUC,CK2VVC,CK2UVC,CK2FCX,  
C  CK2FCY,ZBRWALL,BELV1  
C  MODIFIED TOXIC-ORGANIC CARBON SORPTION OPTIONS ON C44 AND ADDED  
C  ADDITIONAL INPUT FOR POC FRACTIONS ON C45B AND C45C  
C  MOVED SUBGRID SCALE CHANNEL MAPPING TO FOLLOW READ OF MODCHAN.INP  
C  MODIFIED ISDRY OPTION TO INCLUDE TRANSPORT BYPASS  
C  ADDED SED-TOX DEBUG FLAG ISDTXBUG  
C  
      USE GLOBAL  
      USE DRIFTER,ONLY:DRIFTERINP,AREA_CENTRD
#ifdef key_mpi
      USE PARALLEL_MPI
      USE MPI
#endif
      REAL*4       SEEPRATE(1000)
      CHARACTER*80 TEXT,TITLE  
      CHARACTER*10 CDUM  
      CHARACTER*3  NCARD  
      CHARACTER    CCMRM*1, ADUMMY*5                                         !  EJH
      LOGICAL PARSE_LOGICAL, status  
      REAL,ALLOCATABLE,DIMENSION(:)::RMULADS  
      REAL,ALLOCATABLE,DIMENSION(:)::ADDADS
      REAL,ALLOCATABLE,DIMENSION(:)::CONINIT,ATMWHT_TEMP,WNDWHT_TEMP
      REAL,ALLOCATABLE,DIMENSION(:,:)::DXP_GLOB
      REAL,ALLOCATABLE,DIMENSION(:,:)::DYP_GLOB
      REAL,ALLOCATABLE,DIMENSION(:,:)::HMP_GLOB
      INTEGER::ITURB=0,IPMC 
      INTEGER::RI_GLOBAL,RJ_GLOBAL
      INTEGER TINPUT,YBEG,MBEG,DBEG   ! TINPUT is used as YYYYMMDD format
      INTEGER*8:: jd_out
      REAL(8) TBEGIN_IN
      LOGICAL EXIST
      LOGICAL CELL_INSIDE_DOMAIN
      ALLOCATE(RMULADS(NSTM))  
      ALLOCATE(ADDADS(NSTM))  
      ALLOCATE(CONINIT(KCM))
      ALLOCATE(ATMWHT_TEMP(NASERM))
      ALLOCATE(WNDWHT_TEMP(NWSERM))
      ALLOCATE(DXP_GLOB(0:GNX,0:GNY))
      ALLOCATE(DYP_GLOB(0:GNX,0:GNY))
      ALLOCATE(HMP_GLOB(0:GNX,0:GNY))
      RMULADS=0.
      ADDADS=0.
C 
      G=9.81  
      PI=ACOS(-1.0) !3.1415926535898  
      PI2=2.*PI  
    2 FORMAT(A80)  
C  
C **  READ MAIN INPUT FILE EFDC.INP  
C  
      IF (PARTID == MASTER_TASK) THEN
        PRINT *,'READING THE MAIN EFDC CONTROL FILE: EFDC.INP'
      END IF
      OPEN(1,FILE='EFDC.INP',STATUS='OLD', ACTION = 'READ')
C  
C1**  READ TITLE CARD
      NCARD='1'  
      CALL SEEK('C1')  
      READ(1,2) TITLE  
      WRITE(7,1002)NCARD  
      WRITE(7,2) TITLE  
C  
C2**  READ RESTART AND DIAGNOSTIC SWITCHES
      NCARD='2'  
      CALL SEEK('C2')  
      READ(1,*,IOSTAT=ISO) ISRESTI,ISRESTO,ISRESTR,ISPAR,ISLOG,ISDIVEX,  
     &    ISNEGH,ISMMC,ISBAL,IS2TIM,ISHOW  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ISRESTI,ISRESTO,ISRESTR,ISPAR,ISLOG,ISDIVEX,  
     &    ISNEGH,ISMMC,ISBAL,IS2TIM,ISHOW  
      IF(ISMMC.LT.0)THEN
        DEBUG=.TRUE.
        ISMMC=0
      ELSE
        DEBUG=.FALSE.
      ENDIF
      IF(ISO.GT.0) GOTO 100  
C  
C3**  READ RELAXATION PARAMETERS
      NCARD='3'  
      CALL SEEK('C3')  
      READ(1,*,IOSTAT=ISO) RP,RSQM,ITERM,IRVEC,RPADJ,  
     &    RSQMADJ,ITRMADJ,ITERHPM,IDRYCK,ISDSOLV,FILT3TL 
      IF(ITRMADJ.LT.1)ITRMADJ=1 
      WRITE(7,1002)NCARD  
      WRITE(7,*) RP,RSQM,ITERM,IRVEC,RPADJ,  
     &    RSQMADJ,ITRMADJ,ITERHPM,IDRYCK,ISDSOLV,FILT3TL  
      IF(ISO.GT.0) GOTO 100  
      IF(IRVEC.NE. 0.AND.IRVEC.NE.   9.AND.
     &   IRVEC.NE.99.AND.IRVEC.NE.9999)STOP 'INVALID IRVEC'
C     
C4**  READ LONGTERM MASS TRANSPORT INTEGRATION ONLY SWITCHES
      NCARD='4'  
      CALL SEEK('C4')  
      READ(1,*,IOSTAT=ISO) ISLTMT,ISSSMMT,ISLTMTS,ISIA,RPIA,RSQMIA,  
     &    ITRMIA,ISAVEC  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ISLTMT,ISSSMMT,ISLTMTS,ISIA,RPIA,RSQMIA,  
     &    ITRMIA,ISAVEC  
      IF(ISO.GT.0) GOTO 100  
C  
C5**  READ MOMENTUM ADVECTION AND DIFFUSION SWITCHES AND MISC 
      NCARD='5'  
      CALL SEEK('C5')  
      READ(1,*,IOSTAT=ISO) ISCDMA,ISHDMF,ISDISP,ISWASP,ISDRY,  
     &    ISQQ,ISRLID,ISVEG,ISVEGL,ISITB,ISEVER,IINTPG  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ISCDMA,ISHDMF,ISDISP,ISWASP,ISDRY,  
     &    ISQQ,ISRLID,ISVEG,ISVEGL,ISITB,ISEVER,IINTPG  
      IF(ISO.GT.0) GOTO 100  
      IDRYTBP=0  
      IF(ISDRY.LT.0)THEN  
        ISDRY=ABS(ISDRY)  
        IDRYTBP=1  
      ENDIF  
      IF(ISWASP.EQ.99)ISICM=1  
      IF(ISRLID.EQ.1) ISDRY=-1  
      IF(ISWASP.EQ.10)ISRCA=1  
      JSWAVE=0  
C PMC      IS1DCHAN=0  
C PMC  IF(ISCDMA.EQ.10) IS1DCHAN=1  
      ISCOSMIC=0  
C 
      NCARD='6'  
      CALL SEEK('C6')  
      DO N=0,8  
        READ(1,*,IOSTAT=ISO) ISTRAN(N),ISTOPT(N),ISCDCA(N),ISADAC(N),  
     &      ISFCT(N),ISPLIT(N),ISADAH(N),ISADAV(N),ISCI(N),ISCO(N)  
        IF(ISCDCA(N).GE.4) ISCOSMIC=1  
        WRITE(7,1002)NCARD  
        WRITE(7,*) ISTRAN(N),ISTOPT(N),ISCDCA(N),ISADAC(N),  
     &      ISFCT(N),ISPLIT(N),ISADAH(N),ISADAV(N),ISCI(N),ISCO(N)  
        IF(ISO.GT.0) GOTO 100  
      ENDDO  
C  

C7**  READ TIME-RELATED INTEGER PARAMETERS
      NCARD='7' 
      CALL SEEK('C7')  
      READ(1,*,IOSTAT=ISO) NTC,NTSPTC,NLTC,NTTC,NTCPP,NTSTBC,NTCNB,  
     &    NTCVB,NTSMMT,NFLTMT,NDRYSTP  
      WRITE(7,1002)NCARD  
      WRITE(7,*) NTC,NTSPTC,NLTC,NTTC,NTCPP,NTSTBC,NTCNB,  
     &    NTCVB,NTSMMT,NFLTMT,NDRYSTP  
      IF(ISO.GT.0) GOTO 100  
C  
C8**  READ TIME-RELATED REAL PARAMETERS
      NCARD='8'  
      CALL SEEK('C8')  
      READ(1,*,IOSTAT=ISO) TCON,TBEGIN_IN,TIDALP,CF,ISCORV,ISDCCA,
     &    ISCFL,ISCFLM,DTSSFAC  
      WRITE(7,1002)NCARD  
      WRITE(7,*) TCON,TBEGIN_IN,TIDALP,CF,ISCORV,ISDCCA,
     &    ISCFL,ISCFLM,DTSSFAC  
      IF(ISO.GT.0) GOTO 100  
      IF(DTSSFAC.GT.0.0)THEN  
        ISDYNSTP=1  
      ELSE  
        ISDYNSTP=0  
      ENDIF  
      IF(IS2TIM.EQ.0)ISDYNSTP=0 

      !C Fearghal O'Donncha; allow for TINPUT to be of form Julian day
      !or YYYYMMDD 
      IF (TBEGIN_IN > 100000) THEN ! arbitrarily assume of form YYYYMMDD
        TINPUT = INT(TBEGIN_IN)
        CALL Conversion(TINPUT, YBEG, MBEG, DBEG)
       ! Then take base day as day relative to beginning of year 
       ! for backwards compatability with previous versions of EFDC
        TBEGIN = JD_OUT(YBEG,MBEG,DBEG) - JD_OUT(YBEG,01,01) ; ! assuming Jan 1st=day 
        YREF = YBEG
        MREF = MBEG
        DREF = DBEG
        IF (PARTID == 0) THEN
          ! Write to log file details on Begin simulationi
          OPEN(607,FILE='CouplingOutput.log',STATUS='UNKNOWN')
          WRITE(607,1111) YBEG,MBEG,DBEG
        END IF
      ELSE
        YREF = 2019 ! default
        TBEGIN = TBEGIN_IN  ! TBEGIN is specified as number of days since start of year
      END IF
1111   FORMAT('SIMULATION BEGINS FOR DATE: ',I0.4,'-',I0.2,'-',I0.2)
   
C9**  READ SPACE RELATED AND SMOOTHING PARAMETERS
      NCARD='9'  
      CALL SEEK('C9')  
      READ(1,*,IOSTAT=ISO) KC,IC_GLOBAL,JC_GLOBAL,LC_GLOBAL,LVC,
     &    ISCLO,NDM,LDM,ISMASK,ISPGNS,NSHMAX,NSBMAX,WSMH,WSMB  
      WRITE(7,1002)NCARD  
      WRITE(7,*) KC,IC_GLOBAL,JC_GLOBAL,LC_GLOBAL,LVC,ISCLO,NDM,LDM,ISMASK,  
     &    ISPGNS,NSHMAX,NSBMAX,WSMH,WSMB  
      IF(ISO.GT.0) GOTO 100  
      IS2LMC=0  
      IF(KC.LT.0) THEN  
        KC=-KC  
        IS2LMC=1  
      ENDIF  
#ifdef key_mpi
      CALL PARENTGRID
      IC = IC_LORP(XID)
      JC = JC_LORP(YID)
#else
      IC = IC_GLOBAL
      JC = JC_GLOBAL
      DO I = 1,IC_GLOBAL
         DO J = 1,JC_GLOBAL
           XLOC(I) = I
           YLOC(J) = J
         END DO
      END DO
#endif
 111    FORMAT(4I7)
        CLOSE(111)  

C **  DOMAIN DECOMPOSITION CHECKS FOR HORIZONTAL LOOPS  
      LCM2T=LC-2  
      IF(NDM.EQ.1) LDM=LCM2T  
      IF(NDM.GE.2) NCHECK=NDM*LDM  
      IF(NDM.GE.2)THEN  
        IF(NCHECK.NE.LCM2T)THEN  
          WRITE(6,6774)  
          STOP  
        ENDIF  
      ENDIF  
 6774 FORMAT(' INCONSISTENT DOMAIN DECOMPOSITION NDM, LDW ON CARD 9')  
C **  ENDDOMAIN DECOMPOSITION CHECKS FOR HORIZONTAL LOOPS  
      IF(KC.GE.2) ISITB=0  
C  
C10*  READ LAYER THICKNESS IN VERTICAL
      NCARD='10'  
      CALL SEEK('C10')  
      DO K=1,KC  
        READ(1,*,IOSTAT=ISO)KDUM,DZC(K)  
        WRITE(7,1002)NCARD  
        WRITE(7,*)KDUM,DZC(K)  
        IF(ISO.GT.0) GOTO 100  
      ENDDO  
C  
C11*  READ GRID, ROUGHNESS, MASKING AND DEPTH PARAMETERS
      NCARD='11'  
      CALL SEEK('C11')  
      READ(1,*,IOSTAT=ISO) DX,DY,DXYCVT,IMDXDY,ZBRADJ,ZBRCVRT,HMIN,  
     &    HADADJ,HCVRT,HDRY,HWET,BELADJ,BELCVRT  
      WRITE(7,1002)NCARD  
      WRITE(7,*) DX,DY,DXYCVT,IMDXDY,ZBRADJ,ZBRCVRT,HMIN,  
     &    HADADJ,HCVRT,HDRY,HWET,BELADJ,BELCVRT  
      IF(ISO.GT.0) GOTO 100  
C  
C11A* READ TWO-LAYER MOMENTUM FLUX AND CURVATURE ACCELERATION 
C     CORRECTION FACTORS
      NCARD='11A'  
      CALL SEEK('C11A')  
      READ(1,*,IOSTAT=ISO) ICK2COR,CK2UUM,CK2VVM,CK2UVM,CK2UUC,  
     &    CK2VVC,CK2UVC,CK2FCX,CK2FCY  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ICK2COR,CK2UUM,CK2VVM,CK2UVM,CK2UUC,  
     &    CK2VVC,CK2UVC,CK2FCX,CK2FCY  
      IF(ISO.GT.0) GOTO 100  
      IF(ICK2COR.GE.1) THEN  
        IS2LMC=ICK2COR  
      END IF  
C  
C11B* READ CORNER CELL BOTTOM STRESS CORRECTION OPTIONS
      NCARD='11B'  
      CALL SEEK('C11B')  
      READ(1,*,IOSTAT=ISO)ISCORTBC,ISCORTBCD,FSCORTBC
      WRITE(7,1002)NCARD
      WRITE(7,*) ISCORTBC,ISCORTBCD,FSCORTBC
      IF(ISO.GT.0) GOTO 100  
C  
C12*  READ TURBULENT DIFFUSION PARAMETERS
      NCARD='12'  
      CALL SEEK('C12')  
      READ(1,*,IOSTAT=ISO) AHO,AHD,AVO,ABO,AVMX,ABMX,VISMUD,AVCON,
     &                     ZBRWALL,ISAVBMX,ISFAVB,ISINWV
      WRITE(7,1002)NCARD
      WRITE(7,*) AHO,AHD,AVO,ABO,AVMX,ABMX,VISMUD,AVCON,ZBRWALL,
     &           ISAVBMX,ISFAVB,ISINWV
      IF(ISO.GT.0) GOTO 100  
C  
C13*  READ TURBULENCE CLOSURE PARAMETERS
      NCARD='13'  
      CALL SEEK('C13')  
      READ(1,*,IOSTAT=ISO) VKC,CTURB,CTURB2B,CTE1,CTE2,CTE3,QQMIN,  
     &    QQLMIN,DMLMIN  
      WRITE(7,1002)NCARD  
      WRITE(7,*) VKC,CTURB,CTURB2B,CTE1,CTE2,CTE3,QQMIN,  
     &    QQLMIN,DMLMIN  
      IF(ISO.GT.0) GOTO 100  
C  
C14*  READ TIDAL & ATMOSPHERIC FORCING, GROUND WATER 
C     AND SUBGRID CHANNEL PARAMETERS
      NCARD='14'  
      CALL SEEK('C14')  
      READ(1,*,IOSTAT=ISO) MTIDE,NWSER,NASER,ISGWIT,ISCHAN,ISWAVE,  
     &    ITIDASM,ISPERC,ISBODYF,ISPNHYDS  
      WRITE(7,1002)NCARD  
      WRITE(7,*) MTIDE,NWSER,NASER,ISGWIT,ISCHAN,ISWAVE,ITIDASM,  
     &    ISPERC,ISBODYF,ISPNHYDS  
      ISWCBL=0  
      ISWVSD=0  
      IF(ISO.GT.0) GOTO 100  
      IF(ISPERC.GT.0)ISGWIT=3   ! *** DSLLC
      
C14A* READ SAND GRAIN ROUGHNESS
      IF(ISWAVE.EQ.3)THEN
        CALL SEEK('C14A')  
        READ(1,*,IOSTAT=ISO) KSW
      ENDIF
C  
      IF(MTIDE.GT.0)THEN  
C15*  READ PERIODIC FORCING (TIDAL) CONSTITUENT SYMBOLS AND PERIODS
        NCARD='15'  
        CALL SEEK('C15')  
        DO M=1,MTIDE  
          READ(1,*,IOSTAT=ISO) SYMBOL(M),TCP(M)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) SYMBOL(M),TCP(M)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
      ENDIF  
C  
C16*  READ SURFACE ELEVATION OR PRESSURE BOUNDARY CONDITION PARAMETERS
      NCARD='16'  
      CALL SEEK('C16')  
      READ(1,*,IOSTAT=ISO) NPBS,NPBW,NPBE,NPBN,NPFOR,NPFORT,  
     &    NPSER,PDGINIT  
      WRITE(7,1002)NCARD  
      WRITE(7,*) NPBS,NPBW,NPBE,NPBN,NPFOR,NPFORT,NPSER,PDGINIT  
      IF(ISO.GT.0) GOTO 100  
      IF(NPFORT.GE.1.AND.DEBUG)THEN  
        OPEN(2,FILE='TIDALBC.OUT')  
        CLOSE(2,STATUS='DELETE')  
        OPEN(2,FILE='TIDALBC.OUT')  
      ENDIF  
C  
      IF(NPFOR.GT.0)THEN  
C17*  READ PERIODIC FORCING (TIDAL) SURFACE ELEVATION OR
        NCARD='17'  
        CALL SEEK('C17')  
        DO NP=1,NPFOR  
          DO M=1,MTIDE  
            IF(NPFORT.EQ.0)THEN  
              READ(1,*,IOSTAT=ISO)NDUM,CDUM,PFAM(NP,M),PFPH(NP,M)  
              WRITE(7,1002)NCARD  
              WRITE(7,*)NDUM,CDUM,PFAM(NP,M),PFPH(NP,M)  
              IF(ISO.GT.0) GOTO 100  
            ELSE  
              READ(1,*,IOSTAT=ISO)NDUM,CDUM,PFAM(NP,M),PFPH(NP,M)  
              RAD=PI2*PFPH(NP,M)/TCP(M)  
              CPFAM0(NP,M)=PFAM(NP,M)*COS(RAD)  
              SPFAM0(NP,M)=PFAM(NP,M)*SIN(RAD)  
              WRITE(2,2068)NDUM,CDUM,PFAM(NP,M),PFPH(NP,M),  
     &            CPFAM0(NP,M),SPFAM0(NP,M)  
              WRITE(7,1002)NCARD  
              WRITE(7,*)NDUM,CDUM,PFAM(NP,M),PFPH(NP,M),  
     &            CPFAM0(NP,M),SPFAM0(NP,M)  
              IF(ISO.GT.0) GOTO 100  
              READ(1,*,IOSTAT=ISO)NDUM,CDUM,PFAM(NP,M),PFPH(NP,M)  
              RAD=PI2*PFPH(NP,M)/TCP(M)  
              CPFAM1(NP,M)=PFAM(NP,M)*COS(RAD)-CPFAM0(NP,M)  
              SPFAM1(NP,M)=PFAM(NP,M)*SIN(RAD)-SPFAM0(NP,M)  
              WRITE(2,2068)NDUM,CDUM,PFAM(NP,M),PFPH(NP,M),  
     &            CPFAM1(NP,M),SPFAM1(NP,M)  
              WRITE(7,1002)NCARD  
              WRITE(7,*)NDUM,CDUM,PFAM(NP,M),PFPH(NP,M),  
     &            CPFAM1(NP,M),SPFAM1(NP,M)  
              CPFAM2(NP,M)=0.0  
              SPFAM2(NP,M)=0.0  
            ENDIF  
C  
          ENDDO  
        ENDDO  
      ENDIF  
C  
      IF(NPBS.GT.0)THEN  
C18*  READ PERIODIC FORCING (TIDAL) ELEVATION BOUNDARY CONDITIONS 
C     ON SOUTH OPEN BOUNDARIES
        NCARD='18'  
        CALL SEEK('C18')  
        IF(NPFORT.EQ.0)THEN  
          DO L=1,NPBS  
           READ(1,*,IOSTAT=ISO)IPBS_GL(L),JPBS_GL(L),ISPBS_GL(L),NPFORS,NPSERS_GL(L)  
            WRITE(7,1002)NCARD  
            WRITE(7,*) IPBS_GL(L),JPBS_GL(L),ISPBS_GL(L),NPFORS,NPSERS_GL(L)  
            IF(ISO.GT.0) GOTO 100  
            DO M=1,MTIDE  
	        IF(NPFORS.EQ.0) EXIT
              RAD=PI2*PFPH(NPFORS,M)/TCP(M)  
              AMP=G*PFAM(NPFORS,M)  
              PCBS_GL(L,M)=AMP*COS(RAD)  
              PSBS_GL(L,M)=AMP*SIN(RAD)  
            ENDDO  
          ENDDO  
        ELSE  
          DO L=1,NPBS  
            READ(1,*,IOSTAT=ISO) IPBS_GL(L),JPBS_GL(L),ISPBS_GL(L),NPFORS,  
     &          NPSERS_GL(L),TPCOORDS
            WRITE(7,1002)NCARD  
            WRITE(7,*) IPBS_GL(L),JPBS_GL(L),ISPBS_GL(L),NPFORS,NPSERS_GL(L),  
     &          TPCOORDS
            IF(ISO.GT.0) GOTO 100  
            DO M=1,MTIDE  
              PCBS_GL(L,M)=CPFAM0(NPFORS,M)+TPCOORDS*CPFAM1(NPFORS,M)  
     &            +TPCOORDS*TPCOORDS*CPFAM2(NPFORS,M)  
              PSBS_GL(L,M)=SPFAM0(NPFORS,M)+TPCOORDS*SPFAM1(NPFORS,M)  
     &            +TPCOORDS*TPCOORDS*SPFAM2(NPFORS,M)  
            WRITE(2,2069)L,SYMBOL(M),PCBS_GL(L,M),PSBS_GL(L,M),IPBS_GL(L),JPBS_GL(L)  
              PCBS_GL(L,M)=G*PCBS_GL(L,M)  
              PSBS_GL(L,M)=G*PSBS_GL(L,M)  
            ENDDO  
          ENDDO  
        ENDIF  
 2068 FORMAT(I4,3X,A2,5X,E14.4,3E14.5,5X,2I5)  
 2069 FORMAT(I4,3X,A2,5X,2E14.4,5X,2I5)  
      ENDIF  
C  
      IF(NPBW.GT.0)THEN  
C19*  READ PERIODIC FORCING (TIDAL) ELEVATION BOUNDARY CONDITIONS
C     ON WEST OPEN BOUNDARIES
        NCARD='19'  
        CALL SEEK('C19')  
        IF(NPFORT.EQ.0)THEN  
          DO L=1,NPBW  
           READ(1,*,IOSTAT=ISO)IPBW_GL(L),JPBW_GL(L),ISPBW_GL(L),NPFORW,NPSERW_GL(L)  
            WRITE(7,1002)NCARD  
            WRITE(7,*) IPBW_GL(L),JPBW_GL(L),ISPBW_GL(L),NPFORW,NPSERW_GL(L)  
            IF(ISO.GT.0) GOTO 100  
            DO M=1,MTIDE  
	        IF(NPFORW.EQ.0) EXIT
              RAD=PI2*PFPH(NPFORW,M)/TCP(M)  
              AMP=G*PFAM(NPFORW,M)  
              PCBW_GL(L,M)=AMP*COS(RAD)  
              PSBW_GL(L,M)=AMP*SIN(RAD)  
            ENDDO  
          ENDDO  
        ELSE  
          DO L=1,NPBW  
            READ(1,*,IOSTAT=ISO) IPBW_GL(L),JPBW_GL(L),ISPBW_GL(L),NPFORW,  
     &          NPSERW_GL(L),TPCOORDW  
            WRITE(7,1002)NCARD  
            WRITE(7,*) IPBW_GL(L),JPBW_GL(L),ISPBW_GL(L),NPFORW,NPSERW_GL(L),  
     &          TPCOORDW  
            IF(ISO.GT.0) GOTO 100  
            DO M=1,MTIDE  
              PCBW_GL(L,M)=CPFAM0(NPFORW,M)+TPCOORDW*CPFAM1(NPFORW,M)  
     &            +TPCOORDW*TPCOORDW*CPFAM2(NPFORW,M)  
              PSBW_GL(L,M)=SPFAM0(NPFORW,M)+TPCOORDW*SPFAM1(NPFORW,M)  
     &            +TPCOORDW*TPCOORDW*SPFAM2(NPFORW,M)  
            WRITE(2,2069)L,SYMBOL(M),PCBW_GL(L,M),PSBW_GL(L,M),IPBW_GL(L),JPBW_GL(L)  
              PCBW_GL(L,M)=G*PCBW_GL(L,M)  
              PSBW_GL(L,M)=G*PSBW_GL(L,M)  
            ENDDO  
          ENDDO  
        ENDIF  
      ENDIF  
C  
      IF(NPBE.GT.0)THEN  
C20*  READ PERIODIC FORCING (TIDAL)ELEVATION BOUNDARY CONDITIONS
C     ON EAST OPEN BOUNDARIES
        NCARD='20'  
        CALL SEEK('C20')  
        IF(NPFORT.EQ.0)THEN  
          DO L=1,NPBE  
           READ(1,*,IOSTAT=ISO)IPBE_GL(L),JPBE_GL(L),ISPBE_GL(L),NPFORE,NPSERE_GL(L)  
            WRITE(7,1002)NCARD  
            WRITE(7,*) IPBE_GL(L),JPBE_GL(L),ISPBE_GL(L),NPFORE,NPSERE_GL(L)  
            IF(ISO.GT.0) GOTO 100  
            DO M=1,MTIDE  
	        IF(NPFORE.EQ.0) EXIT
              RAD=PI2*PFPH(NPFORE,M)/TCP(M)  
              AMP=G*PFAM(NPFORE,M)  
              PCBE_GL(L,M)=AMP*COS(RAD)  
              PSBE_GL(L,M)=AMP*SIN(RAD)  
            ENDDO  
          ENDDO  
        ELSE  
          DO L=1,NPBE  
            READ(1,*,IOSTAT=ISO) IPBE_GL(L),JPBE_GL(L),ISPBE_GL(L),NPFORE,  
     &          NPSERE_GL(L),TPCOORDE  
            WRITE(7,1002)NCARD  
            WRITE(7,*) IPBE_GL(L),JPBE_GL(L),ISPBE_GL(L),NPFORE,NPSERE_GL(L),  
     &          TPCOORDE
            IF(ISO.GT.0) GOTO 100  
            DO M=1,MTIDE  
              PCBE_GL(L,M)=CPFAM0(NPFORE,M)+TPCOORDE*CPFAM1(NPFORE,M)  
     &            +TPCOORDE*TPCOORDE*CPFAM2(NPFORE,M)  
              PSBE_GL(L,M)=SPFAM0(NPFORE,M)+TPCOORDE*SPFAM1(NPFORE,M)  
     &            +TPCOORDE*TPCOORDE*SPFAM2(NPFORE,M)  
            WRITE(2,2069)L,SYMBOL(M),PCBE_GL(L,M),PSBE_GL(L,M),IPBE_GL(L),JPBE_GL(L)  
              PCBE_GL(L,M)=G*PCBE_GL(L,M)  
              PSBE_GL(L,M)=G*PSBE_GL(L,M)  
            ENDDO  
          ENDDO  
        ENDIF  
      ENDIF  
C  
      IF(NPBN.GT.0)THEN  
C21*  READ PERIODIC FORCING (TIDAL) ELEVATION BOUNDARY CONDITIONS
C     ON NORTH OPEN BOUNDARIES
        NCARD='21'  
        CALL SEEK('C21')  
        IF(NPFORT.EQ.0)THEN  
          DO L=1,NPBN  
           READ(1,*,IOSTAT=ISO)IPBN_GL(L),JPBN_GL(L),ISPBN_GL(L),NPFORN,NPSERN_GL(L)  
            WRITE(7,1002)NCARD  
            WRITE(7,*) IPBN_GL(L),JPBN_GL(L),ISPBN_GL(L),NPFORN,NPSERN_GL(L)  
            IF(ISO.GT.0) GOTO 100  
            DO M=1,MTIDE  
	        IF(NPFORN.EQ.0) EXIT
              RAD=PI2*PFPH(NPFORN,M)/TCP(M)  
              AMP=G*PFAM(NPFORN,M)  
              PCBN_GL(L,M)=AMP*COS(RAD)  
              PSBN_GL(L,M)=AMP*SIN(RAD)  
            ENDDO  
          ENDDO  
        ELSE  
          DO L=1,NPBN  
            READ(1,*,IOSTAT=ISO) IPBN_GL(L),JPBN_GL(L),ISPBN_GL(L),NPFORN,  
     &          NPSERN_GL(L),TPCOORDN
            WRITE(7,1002)NCARD  
            WRITE(7,*) IPBN_GL(L),JPBN_GL(L),ISPBN_GL(L),NPFORN,NPSERN_GL(L),  
     &          TPCOORDN
            IF(ISO.GT.0) GOTO 100  
            DO M=1,MTIDE  
              PCBN_GL(L,M)=CPFAM0(NPFORN,M)+TPCOORDN*CPFAM1(NPFORN,M)  
     &            +TPCOORDN*TPCOORDN*CPFAM2(NPFORN,M)  
              PSBN_GL(L,M)=SPFAM0(NPFORN,M)+TPCOORDN*SPFAM1(NPFORN,M)  
     &            +TPCOORDN*TPCOORDN*SPFAM2(NPFORN,M)  
            WRITE(2,2069)L,SYMBOL(M),PCBN_GL(L,M),PSBN_GL(L,M),IPBN_GL(L),JPBN_GL(L)  
              PCBN_GL(L,M)=G*PCBN_GL(L,M)  
              PSBN_GL(L,M)=G*PSBN_GL(L,M)  
            ENDDO  
          ENDDO  
        ENDIF  
      ENDIF  
C  
      IF(NPFORT.GE.1)THEN  
        CLOSE(2)  
      ENDIF  
C
C21A* READ WATER SURFACE ELEVATION AND VELOCITY DATA ASSIMILATION
C 
      IF(ISWSEDA>0.OR.ISUVDA>0)THEN
        NCARD='21A'
        CALL SEEK('C21A') 
        READ(1,*,IOSTAT=ISO)ISWSEDA,NLWSEDA,ISUVDA,NLUVDA,NUVSER,KCUVDA,KCUVPOW
        WRITE(7,1002)NCARD
        WRITE(7,*) ISWSEDA,NLWSEDA,ISUVDA,NLUVDA,NUVSER,KCUVDA,KCUVPOW
        IF(ISO.GT.0) GOTO 100
C
        DO N=1,NUVSER
	    MUVTLAST(N)=1
        ENDDO
      ENDIF
C
C21B* WATER SURFACE ELEVATION DATA ASSIMILATION
C 
      IF(ISWSEDA>0)THEN
        NCARD='21B'
        CALL SEEK('C21B')
        DO N=1,NLWSEDA                                            
          READ(1,*,IOSTAT=ISO)ICWSEDA_GL(N),JCWSEDA_GL(N),NWSESERA_GL(N),TSWSEDA_GL(N)
        ENDDO
        WRITE(7,1002)NCARD
	  DO N=1,NLWSEDA                                            
          WRITE(7,*)ICWSEDA_GL(N),JCWSEDA_GL(N),NWSESERA_GL(N),TSWSEDA_GL(N)
        ENDDO 
        IF(ISO.GT.0) GOTO 100
      ENDIF
C
C21C* VELOCITY DATA ASSIMILATION
CGR  4/7/10  CREATED SEPERATE ASSIMILATION FACTORS FOR U AND V COMPONENTS, TSUUDA AND TSVVDA 
      IF(ISUVDA.GT.0.0)THEN
        NCARD='21C'
        CALL SEEK('C21C')
        NSPNG=10 !Hard wired for 10 sponge layers. Never an input from Makai.
	  IF(NLUVDA.GT.0.AND.NUVSER.GT.0)THEN
	    DO N=1,NLUVDA                                            
            READ(1,*,IOSTAT=ISO)ICUVDA_GL(N),JCUVDA_GL(N),NUVSERA_GL(N),
     &        TSUUDA_GL(N),TSVVDA_GL(N),NORMDIR_GL(N),FSUVDA_GL(N),IWUVDA_GL(N),
     &        IRVUDA_GL(N),RRUVDA_GL(N)
          ENDDO
          WRITE(7,1002)NCARD
	    DO N=1,NLUVDA                                            
            WRITE(7,*)ICUVDA_GL(N),JCUVDA_GL(N),NUVSERA_GL(N),
     &        TSUUDA_GL(N),TSVVDA_GL(N),NORMDIR_GL(N),FSUVDA_GL(N),IWUVDA_GL(N),
     &        IRVUDA_GL(N),RRUVDA_GL(N)
          ENDDO
          IF(ISO.GT.0) GOTO 100
	  ELSEIF(NLUVDA.GT.0.AND.NUVSER.EQ.0)THEN
          ALLOCATE(TSUVDA(1),TSUVDA_GL(1))
	    N=1                                            
          READ(1,*,IOSTAT=ISO)IDUM,JDUM,NUVSERA_GL(N),
     &        TSUVDA_GL(N),FSUVDA_GL(N),IWUVDA_GL(N),IRVUDA_GL(N),RRUVDA_GL(N)
          NUVSERA(N)=NUVSERA_GL(N)
          WRITE(7,1002)NCARD
          WRITE(7,*)IDUM,JDUM,NUVSERA_GL(N),
     &        TSUVDA_GL(N),FSUVDA_GL(N),IWUVDA_GL(N),IRVUDA_GL(N),RRUVDA_GL(N)
          IF(ISO.GT.0) GOTO 100
	  ENDIF
      ENDIF
C  
C22*  READ NUM OF SEDIMENT AMD TOXICS AND NUM OF CONCENTRATION TIME SERIES
      NCARD='22'  
      CALL SEEK('C22')  
      READ(1,*,IOSTAT=ISO) NTOX,NSED,NSND,NCSER(1),NCSER(2),NCSER(3),  
     &    NCSER(4),NTOXSER,NSEDSER,NSNDSER,ISSBAL  
      WRITE(7,1002)NCARD  
      WRITE(7,*) NTOX,NSED,NSND,NCSER(1),NCSER(2),NCSER(3),  
     &    NCSER(4),NTOXSER,NSEDSER,NSNDSER,ISSBAL  
      IF(ISO.GT.0) GOTO 100  
      MTMP=4  
      DO N=1,NTOX  
        MTMP=MTMP+1  
        MSVTOX(N)=MTMP  
      ENDDO  
      DO N=1,NSED  
        MTMP=MTMP+1  
        MSVSED(N)=MTMP  
      ENDDO  
      DO N=1,NSND  
        MTMP=MTMP+1  
        MSVSND(N)=MTMP  
      ENDDO  
      DO N=1,NTOX  
        M=MSVTOX(N)  
        NCSER(M)=NTOXSER  
      ENDDO  
      DO N=1,NSED  
        M=MSVSED(N)  
        NCSER(M)=NSEDSER  
      ENDDO  
      DO N=1,NSND  
        M=MSVSND(N)  
        NCSER(M)=NSNDSER  
      ENDDO  
      IF(ISTRAN(6).EQ.0.AND.ISTRAN(7).EQ.0) THEN  
        ISSBAL=0  ! *** PMC SINGLE LINE
      ENDIF  
C  
C23*  READ VELOCITY, VOL SOUR/SINK, FLOW CONTROL, & WITHDRAW/RETURN DATA
      NCARD='23'  
      CALL SEEK('C23')  
      READ(1,*,IOSTAT=ISO) NVBS,NUBW,NUBE,NVBN,NQSIJ,NQJPIJ,NQSER,NQCTL,  
     &    NQCTLT,NQWR,NQWRSR,ISDIQ  
      WRITE(7,1002)NCARD  
      WRITE(7,*) NVBS,NUBW,NUBE,NVBN,NQSIJ,NQJPIJ,NQSER,NQCTL,  
     &    NQCTLT,NQWR,NQWRSR,ISDIQ  
      IF(ISO.GT.0) GOTO 100  
C  
      IF(NQSIJ.GT.0)THEN  
C24*  READ VOLUME SOURCE/SINK LOCATIONS, MAGNITUDES, & VOL & CONC SERIES
        NCARD='24'  
        CALL SEEK('C24')  
        DO L=1,NQSIJ  
          READ(1,*,IOSTAT=ISO)IQS_GL(L),JQS_GL(L),QSSE_GL(L),NQSMUL_GL(L),  
     &        NQSMF_GL(L),NQSERQ_GL(L),(NCSERQ_GL(L,K),K=1,4),  
     &        NTOXSRQ_GL(L),NSEDSRQ_GL(L),NSNDSRQ_GL(L),QFACTOR_GL(L)  
          WRITE(7,1002)NCARD  
          WRITE(7,*)IQS_GL(L),JQS_GL(L),QSSE_GL(L),NQSMUL_GL(L),  
     &        NQSMF_GL(L),NQSERQ_GL(L),(NCSERQ_GL(L,K),K=1,4),  
     &        NTOXSRQ_GL(L),NSEDSRQ_GL(L),NSNDSRQ_GL(L),QFACTOR_GL(L)
          IF(ISO.GT.0) GOTO 100  
!          DO K=1,KC  
!            QSS(K,L)=QSSE*DZC(K)  
!          ENDDO  
!          DO N=1,NTOX  
!            M=MSVTOX(N)  
!            NCSERQ(L,M)=NTOXSRQ  
!          ENDDO  
!          DO N=1,NSED  
!            M=MSVSED(N)  
!            NCSERQ(L,M)=NSEDSRQ  
!          ENDDO  
!          DO N=1,NSND  
!            M=MSVSND(N)  
!            NCSERQ(L,M)=NSNDSRQ  
!          ENDDO  
        ENDDO  
C  
C25*  READ TIME CONSTANT VOLUMETRIC SOURCE INFLOW CONCENTRATIONS
C     SAL,TEM,DYE,SFL,TOX(1 TO NOTX)
        NCARD='25'  
        CALL SEEK('C25')  
        MMAX=4+NTOX  
        DO L=1,NQSIJ  
          READ(1,*,IOSTAT=ISO) (CQSE_GL(L,M),M=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CQSE_GL(L,M),M=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
!          DO MS=1,MMAX  
!            DO K=1,KC  
!              CQS(K,L,MS)=CQSE(MS)  
!            ENDDO  
!          ENDDO  
        ENDDO  
C  
C26*  READ TIME CONSTANT VOLUMETRIC SOURCE INFLOW CONCENTRATIONS
C     SED(1 TO NSED),SND(1 TO NSND)
        NCARD='26'  
        CALL SEEK('C26')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NQSIJ  
          READ(1,*,IOSTAT=ISO) (CQSE_GL(L,M),M=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CQSE_GL(L,M),M=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
!          DO MS=MMIN,MMAX  
!            DO K=1,KC  
!              CQS(K,L,MS)=CQSE(MS)  
!            ENDDO  
!          ENDDO  
        ENDDO  
      ENDIF  
C  
      IF(NQJPIJ.GT.0)THEN  
C27*  READ JET/PLUME SOURCE LOCATIONS AND PARAMETERS  
        NCARD='27'  
        CALL SEEK('C27')  
        DO L=1,NQJPIJ  
          READ(1,*,IOSTAT=ISO) IDUM,ICALJP(L),IQJP(L),JQJP(L),KQJP(L),  
     &        NPORTJP(L),XJETL(L),YJETL(L),ZJET(L),PHJET(L),THJET(L),  
     &        DJET(L),CFRD(L),DJPER(L)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) IDUM,ICALJP(L),IQJP(L),JQJP(L),KQJP(L),  
     &        NPORTJP(L),XJETL(L),YJETL(L),ZJET(L),PHJET(L),THJET(L),  
     &        DJET(L),CFRD(L),DJPER(L)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C28*  READ JET/PLUME SOURCE LOCATIONS AND PARAMETERS  
        NCARD='28'  
        CALL SEEK('C28')  
        DO L=1,NQJPIJ  
          READ(1,*,IOSTAT=ISO) IDUM,NJEL(L),NJPMX(L),ISENT(L),ISTJP(L),  
     &        NUDJP(L),IOUTJP(L),NZPRJP(L),ISDJP(L),IUPCJP(L),  
     &        JUPCJP(L),KUPCJP(L)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) IDUM,NJEL(L),NJPMX(L),ISENT(L),ISTJP(L),  
     &        NUDJP(L),IOUTJP(L),NZPRJP(L),ISDJP(L),IUPCJP(L),  
     &        JUPCJP(L),KUPCJP(L)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
!        write(*,*) 'READ C29'
C29*  READ ADDITIONAL JET/PLUME PARAMETERS  
        NCARD='29'  
        CALL SEEK('C29')  
        DO L=1,NQJPIJ  
          READ(1,*,IOSTAT=ISO) IDUM,QQCJP(L),NQSERJP(L),NQWRSERJP(L),  
     &        NCSERJP(L,1),NCSERJP(L,2),NCSERJP(L,3),  
     &        NCSERJP(L,4),NTXSRJP,NSDSRJP,NSNSRJP  
          WRITE(7,1002)NCARD  
          WRITE(7,*) IDUM,QQCJP(L),NQSERJP(L),NQWRSERJP(L),  
     &        NCSERJP(L,1),NCSERJP(L,2),NCSERJP(L,3),  
     &        NCSERJP(L,4),NTXSRJP,NSDSRJP,NSNSRJP  
          NUDJPC(L)=1  
          IF(ISO.GT.0) GOTO 100  
          DO N=1,NTOX  
            M=MSVTOX(N)  
            NCSERJP(L,M)=NTXSRJP  
          ENDDO  
          DO N=1,NSED  
            M=MSVSED(N)  
            NCSERJP(L,M)=NSDSRJP  
          ENDDO  
          DO N=1,NSND  
            M=MSVSND(N)  
            NCSERJP(L,M)=NSNSRJP  
          ENDDO  
          IF(ICALJP(L).EQ.2)THEN  
            QWRCJP(L)=QQCJP(L)  
            QQCJP(L)=0.  
          ELSE  
            QWRCJP(L)=0.  
          ENDIF  
        ENDDO  
        IF(NQJPIJ.GT.1)THEN  
          DO L=2,NQJPIJ  
            NUDJP(L)=NUDJP(1)  
          ENDDO  
        ENDIF  
C  
C30*  READ TIME CONSTANT INFLOW CONCENTRATIONS FOR TIME CONSTANT 
C     JET/PLUME SOURCES SAL,TEM,DYE,SFL,TOX(1 TO NOTX)
        NCARD='30'  
        CALL SEEK('C30')  
        MMAX=4+NTOX  
        DO L=1,NQJPIJ  
          READ(1,*,IOSTAT=ISO) (CQSE(M),M=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CQSE(M),M=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
          IF(ICALJP(L).EQ.1)THEN  
            DO MS=1,MMAX  
              CWRCJP(L,MS)=0.0  
              DO K=1,KC  
                CQCJP(K,L,MS)=CQSE(MS)  
              ENDDO  
            ENDDO  
          ELSE  
            DO MS=1,MMAX  
              CWRCJP(L,MS)=CQSE(MS)  
              DO K=1,KC  
                CQCJP(K,L,MS)=0.0  
              ENDDO  
            ENDDO  
          ENDIF  
        ENDDO  
C  
C31*  READ TIME CONSTANT INFLOW CONCENTRATIONS FOR TIME CONSTANT 
C     JET/PLUME SOURCES SED(1 TO NSED),SND(1 TO NSND)
        NCARD='31'  
        CALL SEEK('C31')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NQJPIJ  
          READ(1,*,IOSTAT=ISO) (CQSE(M),M=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CQSE(M),M=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
          IF(ICALJP(L).EQ.1)THEN  
            DO MS=MMIN,MMAX  
              CWRCJP(L,MS)=0.  
              DO K=1,KC  
                CQCJP(K,L,MS)=CQSE(MS)  
              ENDDO  
            ENDDO  
          ELSE  
            DO MS=MMIN,MMAX  
              CWRCJP(L,MS)=CQSE(MS)  
              DO K=1,KC  
                CQCJP(K,L,MS)=0.  
              ENDDO  
            ENDDO  
          ENDIF  
        ENDDO  
      ENDIF  
C  
      IF(NQCTL.GT.0)THEN  
C32*  READ SURF ELEV OR PRESS DEPENDENT FLOW CONTROL STRUCTURE INFO
        NCARD='32'  
        CALL SEEK('C32')  
        DO L=1,NQCTL  
          READ(1,*,IOSTAT=ISO)IU_TEMP,JU_TEMP,ID_TEMP,JD_TEMP,  
     &        NQCTYP(L),NQCTLQ(L),NQCMUL(L),NQCMFU(L),  
     &        NQCMFD(L),BQCMFU(L),BQCMFD(L)  
               IQCTLU(L) = XLOC(IU_TEMP)
               JQCTLU(L) = YLOC(JU_TEMP)
               IQCTLD(L) = XLOC(ID_TEMP)
               JQCTLD(L) = YLOC(JD_TEMP)
          WRITE(7,1002)NCARD  
          WRITE(7,*)IQCTLU(L),JQCTLU(L),IQCTLD(L),JQCTLD(L),  
     &        NQCTYP(L),NQCTLQ(L),NQCMUL(L),NQCMFU(L),  
     &        NQCMFD(L),BQCMFU(L),BQCMFD(L)  
          IF(ISO.GT.0) GOTO 100  
          DO K=1,KC  
            QCTLTO(K,L)=0.  
            QCTLT(K,L)=0.  
          ENDDO  
        ENDDO  
      ENDIF  
C  
      IF(NQWR.GT.0)THEN  
C33*  READ FLOW WITHDRAWAL, HEAT OR MATERIAL ADDITION, FLOW RETURN DATA
        NCARD='33'  
        CALL SEEK('C33')  
        DO L=1,NQWR  
          READ(1,*,IOSTAT=ISO)IQWRU(L),JQWRU(L),KQWRU(L),  
     &        IQWRD(L),JQWRD(L),KQWRD(L),QWR(L),  
     &        NQWRSERQ(L),NQWRMFU(L),NQWRMFD(L),BQWRMFU(L),BQWRMFD(L),  
     &        ANGWRMFD(L)  
          WRITE(7,1002)NCARD  
          WRITE(7,*)IQWRU(L),JQWRU(L),KQWRU(L),  
     &        IQWRD(L),JQWRD(L),KQWRD(L),QWR(L),  
     &        NQWRSERQ(L),NQWRMFU(L),NQWRMFD(L),BQWRMFU(L),BQWRMFD(L),  
     &        ANGWRMFD(L)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C34*  READ TIME CONSTANT WITHDRAWAL,ADD,RETURN, CONCENTRATION INCREASES
C     SAL,TEM,DYE,SFL,TOX(1 TO NTOX)
        NCARD='34'  
        CALL SEEK('C34')  
        MMAX=4+NTOX  
        DO L=1,NQWR  
          READ(1,*,IOSTAT=ISO) (CQWR(L,MS),MS=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CQWR(L,MS),MS=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C35*  READ TIME CONSTANT WITHDRAWAL,ADD,RETURN, CONCENTRATION INCREASES
C     SED(1 TO NSED),SND(1 TO NSND)
        NCARD='35'  
        CALL SEEK('C35')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NQWR  
          READ(1,*,IOSTAT=ISO) (CQWR(L,MS),MS=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CQWR(L,MS),MS=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
      ENDIF  
C  
      IF(NSED.GT.0.OR.NSND.GT.0)THEN  
C36*  SEDIMENT INITIALIZATION AND WATER COLUMN/BED REPRESENTATION OPTIONS 
        NCARD='36'  
        CALL SEEK('C36')  
        READ(1,*,IOSTAT=ISO)ISEDINT,ISEDBINT,ISEDWC,ISMUD,ISNDWC,ISEDVW,  
     &      ISNDVW,KB,ISDTXBUG  
        WRITE(7,1002)NCARD  
        WRITE(7,*)ISEDINT,ISEDBINT,ISEDWC,ISMUD,ISNDWC,ISEDVW,  
     &      ISNDVW,KB,ISDTXBUG  
        IF(ISO.GT.0) GOTO 100  
C
C36A*  SEDIMENT INITIALIZATION AND WATER COLUMN/BED REPRESENTATION OPTIONS 
        NCARD='36A'  
        CALL SEEK('C36A')  
        COEFTSBL=4.0
        READ(1,*,IOSTAT=ISO)ISBEDSTR,ISBSDFUF,COEFTSBL,VISMUDST
        WRITE(7,1002)NCARD
        WRITE(7,*)ISBEDSTR,ISBSDFUF,COEFTSBL,VISMUDST
        IF(ISO.GT.0) GOTO 100
C
C36B*  SEDIMENT INITIALIZATION AND WATER COLUMN/BED REPRESENTATION OPTIONS 
        NCARD='36B'  
        CALL SEEK('C36B')  
        READ(1,*,IOSTAT=ISO)ISEDAL,ISNDAL,IALTYP,IALSTUP,
     &                     ISEDEFF,HBEDAL,COEHEFF,COEHEFF2
        WRITE(7,1002)NCARD
        WRITE(7,*)ISEDAL,ISNDAL,IALTYP,IALSTUP,
     &                     HBEDAL,COEHEFF,COEHEFF2
        IF(ISO.GT.0) GOTO 100
C  
C37*  BED MECHANICAL PROPERTIES PARAMETER SET 1
        NCARD='37'  
        CALL SEEK('C37')  
        READ(1,*,IOSTAT=ISO)ISEDDT,IBMECH,IMORPH,HBEDMAX,BEDPORC,  
     &      SEDMDMX,SEDMDMN,SEDVDRD,SEDVDRM,SEDVRDT  
        WRITE(7,1002)NCARD  
        WRITE(7,*)ISEDDT,IBMECH,IMORPH,HBEDMAX,BEDPORC,  
     &      SEDMDMX,SEDMDMN,SEDVDRD,SEDVDRM,SEDVRDT  
        IF(ISO.GT.0) GOTO 100  
        ISEDDTC=0  
        IF(IBMECH.EQ.0) THEN  
          SNDVDRD=BEDPORC/(1.-BEDPORC)  
          SEDVDRD=BEDPORC/(1.-BEDPORC)  
          SEDVDRM=SEDVDRD  
        END IF  
        IF(IBMECH.GE.1) THEN  
          ISEDBINT=1  
          IF(ISEDINT.EQ.1) ISEDINT=3  
          ISEDINT=MAX(ISEDINT,2)  
        END IF  
        SNDVDRD=BEDPORC/(1.-BEDPORC)  
        DO NS=1,NSED  
          VDRDEPO(NS)=SEDVDRD  
        ENDDO  
        DO NS=1,NSND  
          NX=NS+NSED  
          VDRDEPO(NX)=SNDVDRD  
        ENDDO  
C  
C38*  BED MECHANICAL PROPERTIES PARAMETER SET 2
        NCARD='38'  
        CALL SEEK('C38')  
        READ(1,*,IOSTAT=ISO)IBMECHK,BMECH1,BMECH2,BMECH3,BMECH4,BMECH5,  
     &      BMECH6  
        WRITE(7,1002)NCARD  
        WRITE(7,*)IBMECHK,BMECH1,BMECH2,BMECH3,BMECH4,BMECH5,BMECH6  
        IF(ISO.GT.0) GOTO 100  
      ENDIF  
C  
      IF(NSED.GT.0)THEN  
C39*  READ COHESIVE SEDIMENT PARAMETER SET 1 REPEAT DATA LINE NSED TIMES
        NCARD='39'  
        CALL SEEK('C39')
        HADJ=0.0  
        DO N=1,NSED  
          READ(1,*,IOSTAT=ISO)SEDO(N),SEDBO(N),SDEN(N),SSG(N),  
     &        WSEDO(N),SEDN(N),SEXP(N),TAUD(N),ISEDSCOR(N)  
          WRITE(7,1002)NCARD  
          WRITE(7,*)SEDO(N),SEDBO(N),SDEN(N),SSG(N),  
     &        WSEDO(N),SEDN(N),SEXP(N),TAUD(N)  
          IF(ISO.GT.0) GOTO 100  
          SEDDIA(N)=0.
          HADJ=SEDN(1)  
        ENDDO  
        IF(HADJ.LT.HWET)HADJ=HWET  ! *** PMC-PROVIDE MORE CONTROL FOR MORPH CHANGE LIMITS
C  
C40*  READ COHESIVE SEDIMENT PARAMETER SET 2 REPEAT DATA LINE NSED TIMES
        NCARD='40'  
        CALL SEEK('C40')  
        DO N=1,NSED  
        READ(1,*,IOSTAT=ISO)IWRSP(N),IWRSPB(N),WRSPO(N),TAUR(N),TAUN(N),  
     &        TEXP(N),VDRRSPO(N),COSEDHID(N)  
          WRITE(7,1002)NCARD  
          WRITE(7,*)IWRSP(N),IWRSPB(N),WRSPO(N),TAUR(N),TAUN(N),TEXP(N),  
     &        VDRRSPO(N),COSEDHID(N)  
          IF(ISO.GT.0) GOTO 100  
C  
          IF(N.EQ.1.AND.IWRSP(N).EQ.999) THEN  
            PRINT *,'READING TAU_CRIT_COH.INP'
            OPEN(1001,FILE='TAU_CRIT_COH.INP',STATUS='OLD')  
            DO L = 2, 4393  
              READ(1001,*,IOSTAT=ISO) (TAUCRCOH(L,K),K=1,10)  
            ENDDO  
            CLOSE(1001)  
          ENDIF  
          IF(ISO.GT.0) GOTO 100  
          ISNDEQ(N)=0  

          ! *** PMC - Mass Erosion is not enabled in EFDC at this time, so ensure disabled
          IF(IWRSPB(N).GT.0)THEN
            PRINT *,' *** WARNING: COHESIVE MASS/BULK EROSION IS NOT ENA
     &BLED IN EFDC!'
            IWRSPB(N)=0
          ENDIF
        ENDDO  
      ENDIF  
C 
      IF(NSND.GT.0)THEN  
C41*  READ NONCOHESIVE SEDIMENT PARAMETER SET 1 REPEAT DATA LINE NSND TIMES  
        NCARD='41'  
        CALL SEEK('C41')  
        DO NX=1,NSND  
          N=NX+NSED  
          READ(1,*,IOSTAT=ISO)SEDO(N),SEDBO(N),SDEN(N),SSG(N),SEDDIA(N),  
     &        WSEDO(N),SEDN(N),SEXP(N),TAUD(N),ISEDSCOR(N)  
          ! *** IF SETTLING VELOCITY IS NEGATIVE, COMPUTE USING VAN RIJN'S FORMULA  
          IF(WSEDO(N).LT.0.0)THEN  
            WSEDO(N)=SETSTVEL(SEDDIA(N),SSG(N))  
          ENDIF  
          WRITE(7,1002)NCARD  
          WRITE(7,*)SEDO(N),SEDBO(N),SDEN(N),SSG(N),SEDDIA(N),  
     &        WSEDO(N),SEDN(N),SEXP(N),TAUD(N),ISEDSCOR(N)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C42*  READ NONCOHESIVE SEDIMENT PARAMETER SET 2 REPEAT DATA LINE NSND TIMES  
        NCARD='42'  
        CALL SEEK('C42')  
        DO NX=1,NSND  
          N=NX+NSED  
          READ(1,*,IOSTAT=ISO)ISNDEQ(N),ISBDLD(N),TAUR(N),TAUN(N),
     &                     TCSHIELDS(N),ISLTAUC(N),IBLTAUC(N),
     &                     IROUSE(NX),ISNDM1(NX),ISNDM2(NX),RSNDM(NX)
C  
C IF TAUR(N) IS NEGATIVE, COMPUTE USING VAN RIJN'S FORMULA  
C       TAUR:     CRITICAL SHIELDS STRESS IN (M/S)**2   (ISNDEQ=2)  
C       TAUN:     EQUAL TO TAUR FOR NONCHOESIVE SED TRANS  (ISNDEQ=2)  
C       TEXP:     CRITICAL SHIELDS PARAMETER  (ISNDEQ=2)  
C  
          DSTR=0.0  
          USTR=0.0  
          IF(TAUR(N).LT.0.0)THEN  
           CALL SETSHLD(TAUR(N),TCSHIELDS(N),SEDDIA(N),SSG(N),DSTR,USTR)  
            TAUN(N)=TAUR(N)  
          ENDIF  
C  
C IF TAUR(N) IS NEGATIVE, COMPUTE USING VAN RIJN'S FORMULA  
C  
          WRITE(7,1002)NCARD  
          WRITE(7,*)ISNDEQ(N),TAUR(N),TAUN(N),TCSHIELDS(N),SEDDIA(N),  
     &        SSG(N),DSTR,USTR  
          IF(ISO.GT.0) GOTO 100  
          IWRSP(N)=0  
          WRSPO(N)=0  
        ENDDO  
C  
C42A*  READ NONCOHESIVE SEDIMENT BED LOAD PARAMETERS 
        NCARD='42A'  
        CALL SEEK('C42A')  
        READ(1,*,IOSTAT=ISO)ISBDLDBC,SBDLDA,SBDLDB,SBDLDG1,  
     &      SBDLDG2,SBDLDG3,SBDLDG4,SBDLDP,ISBLFUC,BLBSNT  
        WRITE(7,1002)NCARD  
        WRITE(7,*)ISBDLDBC,SBDLDA,SBDLDB,SBDLDG1,SBDLDG2,SBDLDG3,  
     &      SBDLDG4,SBDLDP,ISBLFUC,BLBSNT  
        IF(ISO.GT.0) GOTO 100  
      ENDIF  
C  
      IF(NTOX.GT.0)THEN  
C43*  READ TOXIC CONTAMINANT INITIAL CONDITIONS AND PARAMETERS
        NCARD='43'  
        CALL SEEK('C43')  
	  IF(IWRSP(1).LT.90)THEN  ! *** SEDZLJ Edit
          DO NT=1,NTOX  
            READ(1,*,IOSTAT=ISO)NDUM,ITXINT(NT),ITXBDUT(NT),TOXINTW(NT),  
     &          TOXINTB(NT),RKTOXW(NT),TKTOXW(NT),RKTOXB(NT),TRTOXB(NT)  
            WRITE(7,1002)NCARD  
            WRITE(7,*)NDUM,ITXINT(NT),ITXBDUT(NT),TOXINTW(NT),  
     &          TOXINTB(NT),RKTOXW(NT),TKTOXW(NT),RKTOXB(NT),TRTOXB(NT)  
            IF(ISO.GT.0) GOTO 100  
          ENDDO  
        ENDIF
C  
C44*  READ TOXIC CONTAMINANT PARAMETERS
        NCARD='44'  
        CALL SEEK('C44')  
	  IF(IWRSP(1).LT.90)THEN  ! *** SEDZLJ Edit
          DO NT=1,NTOX  
            READ(1,*,IOSTAT=ISO)NDUM,ISTOC(NT),VOLTOX(NT),RMOLTX(NT),  
     &        RKTOXP(NT),SKTOXP(NT),DIFTOX(NT),  
     &        DIFTOXS(NT),PDIFTOX (NT),DPDIFTOX(NT)  
            WRITE(7,1002)NCARD  
            WRITE(7,*)NDUM,ISTOC(NT),VOLTOX(NT),RMOLTX(NT),RKTOXP(NT),  
     &        SKTOXP(NT),DIFTOX(NT),  
     &        DIFTOXS(NT),PDIFTOX (NT),DPDIFTOX(NT)  
            IF(ISO.GT.0) GOTO 100  
	      ISDIFBW(NT)=0
	      IF(DIFTOXS(NT).LT.0.0)THEN
	      DIFTOXS(NT)=ABS(DIFTOXS(NT))
	      ISDIFBW(NT)=1
	      ENDIF
          ENDDO  
        ENDIF
C  
C45*  READ TOXIC CONTAMINANT-SEDIMENT INTERACTION PARAMETERS
        NCARD='45'  
        CALL SEEK('C45')  
	  IF(IWRSP(1).LT.90)THEN  ! *** SEDZLJ Edit
        DO NT=1,NTOX  
          IF(NSED.GT.0)THEN  
            DO N=1,NSED  
              READ(1,*,IOSTAT=ISO)NDUM1,NDUM2,  
     &            ITXPARW(N,NT),TOXPARW(N,NT),CONPARW(N,NT),  
     &            ITXPARB(N,NT),TOXPARB(N,NT),CONPARB(N,NT)  
              WRITE(7,1002)NCARD  
              WRITE(7,*)NDUM1,NDUM2,  
     &            ITXPARW(N,NT),TOXPARW(N,NT),CONPARW(N,NT),  
     &            ITXPARB(N,NT),TOXPARB(N,NT),CONPARB(N,NT)  
              IF(ISO.GT.0) GOTO 100  
            ENDDO  
          ENDIF  
          IF(NSND.GT.0)THEN  
            DO NX=1,NSND  
              N=NX+NSED  
              READ(1,*,IOSTAT=ISO)NDUM1,NDUM2,  
     &            ITXPARW(N,NT),TOXPARW(N,NT),CONPARW(N,NT),  
     &            ITXPARB(N,NT),TOXPARB(N,NT),CONPARB(N,NT)  
              WRITE(7,1002)NCARD  
              WRITE(7,*)NDUM1,NDUM2,  
     &            ITXPARW(N,NT),TOXPARW(N,NT),CONPARW(N,NT),  
     &            ITXPARB(N,NT),TOXPARB(N,NT),CONPARB(N,NT)  
              IF(ISO.GT.0) GOTO 100  
            ENDDO  
          ENDIF  
        ENDDO  
      ENDIF
C  
!        write(*,*) 'READ C45'
C45A*  READ TOXIC CONTAMINANT ORGANIC CARBON PARAMETERS
        NCARD='45A'  
        CALL SEEK('C45A')   
	  IF(IWRSP(1).LT.90.AND.NTOX.GT.0)THEN ! SEDZLJ edit
        READ(1,*,IOSTAT=ISO)ISTDOCW,ISTPOCW,ISTDOCB,ISTPOCB,  
     &      STDOCWC,STPOCWC,STDOCBC,STPOCBC  
        WRITE(7,1002)NCARD  
        WRITE(7,*)ISTDOCW,ISTPOCW,ISTDOCB,ISTPOCB,  
     &      STDOCWC,STPOCWC,STDOCBC,STPOCBC  
        IF(ISO.GT.0) GOTO 100
	  ENDIF  
C  
C45B* READ TOXIC CONTAMINANT-ORGANIC CARBON INTERACTION PARAMETERS
        NCARD='45B'  
        CALL SEEK('C45B')  
	  IF(IWRSP(1).LT.90)THEN  ! *** SEDZLJ Edit
          DO NT=1,NTOX  
            READ(1,*,IOSTAT=ISO)NDUM1,NDUM2,  
     &        ITXPARWC(1,NT),TOXPARWC(1,NT),CONPARWC(1,NT),  
     &        ITXPARBC(1,NT),TOXPARBC(1,NT),CONPARBC(1,NT)  
            WRITE(7,1002)NCARD  
            WRITE(7,*)NDUM1,NDUM2,  
     &        ITXPARWC(1,NT),TOXPARWC(1,NT),CONPARWC(1,NT),  
     &        ITXPARBC(1,NT),TOXPARBC(1,NT),CONPARBC(1,NT)  
            IF(ISO.GT.0) GOTO 100  
            READ(1,*,IOSTAT=ISO)NDUM1,NDUM2,  
     &        ITXPARWC(2,NT),TOXPARWC(2,NT),CONPARWC(2,NT),  
     &        ITXPARBC(2,NT),TOXPARBC(2,NT),CONPARBC(2,NT)  
            WRITE(7,1002)NCARD  
            WRITE(7,*)NDUM1,NDUM2,  
     &        ITXPARWC(2,NT),TOXPARWC(2,NT),CONPARWC(2,NT),  
     &        ITXPARBC(2,NT),TOXPARBC(2,NT),CONPARBC(2,NT)  
            IF(ISO.GT.0) GOTO 100  
          ENDDO  
        ENDIF
C  
C45C* READ TOXIC CONTAMINANT-ORGANIC CARBON WATER COLUMN POC FRACTIONS
        NCARD='45C'  
        CALL SEEK('C45C')  
        WRITE(7,1002)NCARD  
        NTMP=NSED+NSND  
	  IF(IWRSP(1).LT.90)THEN  ! *** SEDZLJ Edit
          DO NT=1,NTOX  
            READ(1,*,IOSTAT=ISO)NDUM1,(FPOCWST(NS,NT),NS=1,NTMP)  
            IF(ISO.GT.0) GOTO 100  
            WRITE(7,*)NDUM1,(FPOCWST(NS,NT),NS=1,NTMP)  
          ENDDO  
	  ENDIF
C  
C RESET INORGANIC SEDIMENT PARTITION COEFFICIENTS BASED ON  
C FRACTION OF POC ASSIGNED TO EACH SEDIMENT CLASS IN WATER COLUMN
C  
        DO NT=1,NTOX  
          IF(ISTOC(NT).GE.2)THEN  
            IF(NSED.GT.0)THEN  
              DO NS=1,NSED  
                ITXPARW(NS,NT)=0  
                TOXPARW(NS,NT)=TOXPARWC(2,NT)  
                CONPARW(NS,NT)=0.  
              ENDDO  
            ENDIF  
            IF(NSND.GT.0)THEN  
              DO NX=1,NSND  
                NS=NSED+NX  
                ITXPARW(NS,NT)=0  
                TOXPARW(NS,NT)=TOXPARWC(2,NT)  
                CONPARW(NS,NT)=0.  
              ENDDO  
            ENDIF  
          ENDIF  
        ENDDO  
C  
C45D* READ TOXIC CONTAMINANT-ORGANIC CARBON SED BED POC FRACTIONS
        NCARD='45D'  
        CALL SEEK('C45D')  
        WRITE(7,1002)NCARD  
        NTMP=NSED+NSND  
	  IF(IWRSP(1).LT.90)THEN  ! *** SEDZLJ Edit
          DO NT=1,NTOX  
            READ(1,*,IOSTAT=ISO)NDUM1,(FPOCBST(NS,NT),NS=1,NTMP)  
            IF(ISO.GT.0) GOTO 100  
            WRITE(7,*)NDUM1,(FPOCBST(NS,NT),NS=1,NTMP)  
          ENDDO  
        ENDIF
C  
C RESET INORGANIC SEDIMENT PARTITION COEFFICIENTS BASED ON  
C FRACTION OF POC ASSIGNED TO EACH SEDIMENT CLASS IN SEDIMENT BED  
C  
        DO NT=1,NTOX  
          IF(ISTOC(NT).GE.2)THEN  
            IF(NSED.GT.0)THEN  
              DO NS=1,NSED  
                ITXPARB(NS,NT)=0  
                TOXPARB(NS,NT)=TOXPARBC(2,NT)  
                CONPARB(NS,NT)=0  
              ENDDO  
            ENDIF  
            IF(NSND.GT.0)THEN  
              DO NX=1,NSND  
                NS=NSED+NX  
                ITXPARB(NS,NT)=0  
                TOXPARB(NS,NT)=TOXPARBC(2,NT)  
                CONPARB(NS,NT)=0  
              ENDDO  
            ENDIF  
          ENDIF  
        ENDDO  
      ENDIF  
C  
C46*  READ BUOYANCY, TEMPERATURE, DYE DATA AND CONCENTRATION BC DATA
      NCARD='46'  
      CALL SEEK('C46')  
      READ(1,*,IOSTAT=ISO)BSC,TEMO,HEQT,RKDYE,NCBS,NCBW,NCBE,NCBN  
      WRITE(7,1002)NCARD  
      WRITE(7,*)BSC,TEMO,HEQT,RKDYE,NCBS,NCBW,NCBE,NCBN  
      IF(ISO.GT.0) GOTO 100  
      IF(BSC.EQ.2)THEN  
        BSC=1.  
        IBSC=1  
      ELSE  
        IBSC=0  
      ENDIF  
      IF(TEMO.LT.0.0)THEN
        TEMO=ABS(TEMO)
        INITTEMP=1
      ELSE
        INITTEMP=0
      ENDIF
C  
      IF(NCBS.GT.0)THEN  
C47*  READ LOCATIONS OF CONC BC'S ON SOUTH BOUNDARIES
        NCARD='47'  
        CALL SEEK('C47')  
        DO L=1,NCBS  
          READ(1,*,IOSTAT=ISO) ICBS_GL(L),JCBS_GL(L),NTSCRS_GL(L),  
     &        (NCSERS_GL(L,K),K=1,4),NTOXSRC,NSEDSRC,NSNDSRC  
          WRITE(7,1002)NCARD  
          WRITE(7,*) ICBS_GL(L),JCBS_GL(L),NTSCRS_GL(L),  
     &        (NCSERS_GL(L,K),K=1,4),NTOXSRC,NSEDSRC,NSNDSRC  
          IF(ISO.GT.0) GOTO 100  
          DO N=1,NTOX  
            M=MSVTOX(N)  
            NCSERS_GL(L,M)=NTOXSRC  
          ENDDO  
          DO N=1,NSED  
            M=MSVSED(N)  
            NCSERS_GL(L,M)=NSEDSRC  
          ENDDO  
          DO N=1,NSND  
            M=MSVSND(N)  
            NCSERS_GL(L,M)=NSNDSRC  
          ENDDO  
        ENDDO  
C  
C48*  READ CONSTANT BOTTOM CONCENTRATION ON SOUTH CONC BOUNDARIES
C     SAL,TEM,DYE,SFL,TOX(1 TO NTOX)
        NCARD='48'  
        CALL SEEK('C48')  
        MMAX=4+NTOX  
        DO L=1,NCBS  
          READ(1,*,IOSTAT=ISO) (CBS_GL(L,1,M),M=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBS_GL(L,1,M),M=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C49*  READ CONSTANT BOTTOM CONCENTRATION ON SOUTH CONC BOUNDARIES
C     SED(1 TO NSED),SND(1,NSND)
        NCARD='49'  
        CALL SEEK('C49')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NCBS  
          READ(1,*,IOSTAT=ISO) (CBS_GL(L,1,M),M=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBS_GL(L,1,M),M=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C50*  READ CONSTANT SURFACE CONCENTRATION ON SOUTH CONC BOUNDARIES
C     SAL,TEM,DYE,SFL,TOX(1 TO NTOX)
        NCARD='50'  
        CALL SEEK('C50')  
        MMAX=4+NTOX  
        DO L=1,NCBS  
          READ(1,*,IOSTAT=ISO) (CBS_GL(L,2,M),M=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBS_GL(L,2,M),M=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C51*  READ CONSTANT SURFACE CONCENTRATION ON SOUTH CONC BOUNDARIES
C     SED(1 TO NSED),SND(1,NSND)
        NCARD='51'  
        CALL SEEK('C51')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NCBS  
          READ(1,*,IOSTAT=ISO) (CBS_GL(L,2,M),M=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBS_GL(L,2,M),M=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
      ENDIF  
C  
      IF(NCBW.GT.0)THEN  
C52*  READ LOCATIONS OF CONC BC'S ON WEST BOUNDARIES
        NCARD='52'  
        CALL SEEK('C52')  
        DO L=1,NCBW  
          READ(1,*,IOSTAT=ISO) ICBW_GL(L),JCBW_GL(L),NTSCRW_GL(L),  
     &        (NCSERW_GL(L,K),K=1,4),NTOXSRC,NSEDSRC,NSNDSRC  
          WRITE(7,1002)NCARD  
          WRITE(7,*) ICBW_GL(L),JCBW_GL(L),NTSCRW_GL(L),  
     &        (NCSERW_GL(L,K),K=1,4),NTOXSRC,NSEDSRC,NSNDSRC  
          IF(ISO.GT.0) GOTO 100  
          DO N=1,NTOX  
            M=MSVTOX(N)  
            NCSERW_GL(L,M)=NTOXSRC  
          ENDDO  
          DO N=1,NSED  
            M=MSVSED(N)  
            NCSERW_GL(L,M)=NSEDSRC  
          ENDDO  
          DO N=1,NSND  
            M=MSVSND(N)  
            NCSERW_GL(L,M)=NSNDSRC  
          ENDDO  
        ENDDO  
C  
C53*  READ CONSTANT BOTTOM CONCENTRATION ON WEST CONC BOUNDARIES
C     SAL,TEM,DYE,SFL,TOX(1 TO NTOX)
        NCARD='53'  
        CALL SEEK('C53')  
        MMAX=4+NTOX  
        DO L=1,NCBW  
          READ(1,*,IOSTAT=ISO) (CBW_GL(L,1,M),M=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBW_GL(L,1,M),M=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C54*  READ CONSTANT BOTTOM CONCENTRATION ON WEST CONC BOUNDARIES
C     SED(1 TO NSED),SND(1,NSND)
        NCARD='54'  
        CALL SEEK('C54')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NCBW  
          READ(1,*,IOSTAT=ISO) (CBW_GL(L,1,M),M=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBW_GL(L,1,M),M=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C55*  READ CONSTANT SURFACE CONCENTRATION ON WEST CONC BOUNDARIES
C     SAL,TEM,DYE,SFL,TOX(1 TO NTOX)
        NCARD='55'  
        CALL SEEK('C55')  
        MMAX=4+NTOX  
        DO L=1,NCBW  
          READ(1,*,IOSTAT=ISO) (CBW_GL(L,2,M),M=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBW_GL(L,2,M),M=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C56*  READ CONSTANT SURFACE CONCENTRATION ON WEST CONC BOUNDARIES
C     SED(1 TO NSED),SND(1,NSND)
        NCARD='56'  
        CALL SEEK('C56')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NCBW  
          READ(1,*,IOSTAT=ISO) (CBW_GL(L,2,M),M=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBW_GL(L,2,M),M=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
      ENDIF  
C  
      IF(NCBE.GT.0)THEN  
C57*  READ LOCATIONS OF CONC BC'S ON EAST BOUNDARIES
        NCARD='57'  
        CALL SEEK('C57')  
        DO L=1,NCBE  
          READ(1,*,IOSTAT=ISO) ICBE_GL(L),JCBE_GL(L),NTSCRE_GL(L),  
     &        (NCSERE_GL(L,K),K=1,4),NTOXSRC,NSEDSRC,NSNDSRC  
          WRITE(7,1002)NCARD  
          WRITE(7,*) ICBE_GL(L),JCBE_GL(L),NTSCRE_GL(L),  
     &        (NCSERE(L,K),K=1,4),NTOXSRC,NSEDSRC,NSNDSRC  
          IF(ISO.GT.0) GOTO 100  
          DO N=1,NTOX  
            M=MSVTOX(N)  
            NCSERE_GL(L,M)=NTOXSRC  
          ENDDO  
          DO N=1,NSED  
            M=MSVSED(N)  
            NCSERE_GL(L,M)=NSEDSRC  
          ENDDO  
          DO N=1,NSND  
            M=MSVSND(N)  
            NCSERE_GL(L,M)=NSNDSRC  
          ENDDO  
        ENDDO  
C  
C58*  READ CONSTANT BOTTOM CONCENTRATION ON EAST CONC BOUNDARIES
C     SAL,TEM,DYE,SFL,TOX(1 TO NTOX)
        NCARD='58'  
        CALL SEEK('C58')  
        MMAX=4+NTOX  
        DO L=1,NCBE  
          READ(1,*,IOSTAT=ISO) (CBE_GL(L,1,M),M=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBE_GL(L,1,M),M=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C59*  READ CONSTANT BOTTOM CONCENTRATION ON EAST CONC BOUNDARIES
C     SED(1 TO NSED),SND(1,NSND)
        NCARD='59'  
        CALL SEEK('C59')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NCBE  
          READ(1,*,IOSTAT=ISO) (CBE_GL(L,1,M),M=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBE_GL(L,1,M),M=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C60*  READ CONSTANT SURFACE CONCENTRATION ON EAST CONC BOUNDARIES
C     SAL,TEM,DYE,SFL,TOX(1 TO NTOX)
        NCARD='60'  
        CALL SEEK('C60')  
        MMAX=4+NTOX  
        DO L=1,NCBE  
          READ(1,*,IOSTAT=ISO) (CBE_GL(L,2,M),M=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBE_GL(L,2,M),M=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C61*  READ CONSTANT SURFACE CONCENTRATION ON EAST CONC BOUNDARIES
C     SED(1 TO NSED),SND(1,NSND)
        NCARD='61'  
        CALL SEEK('C61')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NCBE  
          READ(1,*,IOSTAT=ISO) (CBE_GL(L,2,M),M=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBE_GL(L,2,M),M=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
      ENDIF  
C  
      IF(NCBN.GT.0)THEN  
C62*  READ LOCATIONS OF CONC BC'S ON NORTH BOUNDARIES
        NCARD='62'  
        CALL SEEK('C62')  
        DO L=1,NCBN  
          READ(1,*,IOSTAT=ISO) ICBN_GL(L),JCBN_GL(L),NTSCRN_GL(L),  
     &        (NCSERN_GL(L,K),K=1,4),NTOXSRC,NSEDSRC,NSNDSRC  
          WRITE(7,1002)NCARD  
          WRITE(7,*) ICBN_GL(L),JCBN_GL(L),NTSCRN_GL(L),  
     &        (NCSERN(L,K),K=1,4),NTOXSRC,NSEDSRC,NSNDSRC  
          IF(ISO.GT.0) GOTO 100  
          DO N=1,NTOX  
            M=MSVTOX(N)  
            NCSERN_GL(L,M)=NTOXSRC  
          ENDDO  
          DO N=1,NSED  
            M=MSVSED(N)  
            NCSERN_GL(L,M)=NSEDSRC  
          ENDDO  
          DO N=1,NSND  
            M=MSVSND(N)  
            NCSERN_GL(L,M)=NSNDSRC  
          ENDDO  
        ENDDO  
C  
C63*  READ CONSTANT BOTTOM CONCENTRATION ON NORTH CONC BOUNDARIES
C     SAL,TEM,DYE,SFL,TOX(1 TO NTOX)
        NCARD='63'  
        CALL SEEK('C63')  
        MMAX=4+NTOX  
        DO L=1,NCBN  
          READ(1,*,IOSTAT=ISO) (CBN_GL(L,1,M),M=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBN_GL(L,1,M),M=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C64*  READ CONSTANT BOTTOM CONCENTRATION ON NORTH CONC BOUNDARIES
C     SED(1 TO NSED),SND(1,NSND)
        NCARD='64'  
        CALL SEEK('C64')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NCBN  
          READ(1,*,IOSTAT=ISO) (CBN_GL(L,1,M),M=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBN_GL(L,1,M),M=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C65*  READ CONSTANT SURFACE CONCENTRATION ON NORTH CONC BOUNDARIES
C     SAL,TEM,DYE,SFL,TOX(1 TO NTOX)
        NCARD='65'  
        CALL SEEK('C65')  
        MMAX=4+NTOX  
        DO L=1,NCBN  
          READ(1,*,IOSTAT=ISO) (CBN_GL(L,2,M),M=1,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBN_GL(L,2,M),M=1,MMAX)  
          IF(ISO.GT.0) GOTO 100  
      ENDDO  
C  
C66*  READ CONSTANT SURFACE CONCENTRATION ON NORTH CONC BOUNDARIES
C     SED(1 TO NSED),SND(1,NSND)
        NCARD='66'  
        CALL SEEK('C66')  
        MMIN=MMAX+1  
        MMAX=MMAX+NSED+NSND  
        DO L=1,NCBN  
          READ(1,*,IOSTAT=ISO) (CBN_GL(L,2,M),M=MMIN,MMAX)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) (CBN_GL(L,2,M),M=MMIN,MMAX)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
      ENDIF  
C  
C66A*  READ CONCENTRATION DATA ASSIMILATION PARAMETERS
      NCARD='66A'  
      CALL SEEK('C66A')
      READ(1,*,IOSTAT=ISO) NLCDA,TSCDA,(ISCDA(K),K=1,7)  
      WRITE(7,1002)NCARD  
      WRITE(7,*) NLCDA,TSCDA,(ISCDA(K),K=1,7)  
      IF(ISO.GT.0) GOTO 100  
C
C66B*  READ CONCENTRATION DATA ASSIMILATION LOCATIONS AND
C      SERIES IDENTIFIERS
      IF(NLCDA.GT.0.AND.NUVSER==0)THEN  !SCJ check this it is problematic
        NCARD='66B'  
        CALL SEEK('C66B')  
        WRITE(7,1002)NCARD
        DO L=1,NLCDA
         READ(1,*,IOSTAT=ISO) ITPCDA(L),ICDA(L),JCDA(L),
     &                           ICCDA(L),JCCDA(L),(NCSERA(L,K),K=1,7) !Makai added variables KCSERDA,KCSERPOW K=1->salinity, K=2->temperature, K=3->dye, K=4->shellfish, K=5->toxics, K=6->cohesive sediment, K=7->noncohesive sediment
         KCSERDA=0
         KCSERPOW=0
         WRITE(7,*)  ITPCDA(L),ICDA(L),JCDA(L),
     &                           ICCDA(L),JCCDA(L),(NCSERA(L,K),K=1,7)
         IF(ISO.GT.0) GOTO 100
        ENDDO
        PAUSE 'Needs to be updated for MPI'
      ENDIF
C
C67*  READ NEUTRALLY BUOYANT PARTICLE DRIFTER DATA
      NCARD='67'  
      CALL SEEK('C67')  
      READ(1,*,IOSTAT=ISO) ISPD,NPD,NPDRT,NWPD,ISLRPD,ILRPD1,ILRPD2,  
     &    JLRPD1, JLRPD2, MLRPDRT,IPLRPD  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ISPD,NPD,NPDRT,NWPD,ISLRPD,ILRPD1,ILRPD2,  
     &    JLRPD1, JLRPD2, MLRPDRT,IPLRPD  
      IF(ISO.GT.0) GOTO 100  
C  
      IF(NPD.GT.0)THEN  
C68*  READ NEUTRALLY BUOYANT PARTICLE INITIAL POSITIONS
        NCARD='68'  
        CALL SEEK('C68')  
        DO NP=1,NPD  
          READ(1,*,IOSTAT=ISO) RI_GLOBAL,RJ_GLOBAL,RK(NP)
          ! MAP FROM THE GLOBAL DOMAIN FOR COMPUTATION (RK DOES NOT CHANGE)
          I = XLOC(RI_GLOBAL)
          J = YLOC(RJ_GLOBAL)
          IF (I > 0 .AND. I <=IC) THEN
            RI(NP) = I
          END IF
          IF (J > 0 .AND. J <=JC) THEN
            RJ(NP) = J
          END IF
          WRITE(7,1002)NCARD  
          WRITE(7,*)  RI(NP),RJ(NP),RK(NP)
        ENDDO  
      ENDIF  
C  
C69*  CONSTANTS FOR LONGITUDE AND LATITUDE OF CELL CENTERS
      NCARD='69'  
      CALL SEEK('C69')  
      READ(1,*,IOSTAT=ISO) CDLON1,CDLON2,CDLON3,CDLAT1,CDLAT2,CDLAT3  
      WRITE(7,1002)NCARD  
      WRITE(7,*) CDLON1,CDLON2,CDLON3,CDLAT1,CDLAT2,CDLAT3  
      IF(ISO.GT.0) GOTO 100  
C  
C70*  CONTROLS FOR WRITING ASCII OR BINARY DUMP FILES 
      NCARD='70'  
      CALL SEEK('C70')  
      READ(1,*,IOSTAT=ISO)ISDUMP,ISADMP,NSDUMP,TSDUMP,TEDUMP,ISDMPP,  
     &    ISDMPU,ISDMPW,ISDMPT,IADJDMP  
      WRITE(7,1002)NCARD  
      WRITE(7,*)ISDUMP,ISADMP,NSDUMP,TSDUMP,TEDUMP,ISDMPP,  
     &    ISDMPU,ISDMPW,ISDMPT,IADJDMP  
      IF(ISO.GT.0) GOTO 100  
      JSDUMP=1  
      NCDUMP=1  
C  
C71*  CONTROLS FOR HORIZONTAL PLANE SCALAR FIELD CONTOURING 
      NCARD='71'  
      CALL SEEK('C71')  
      DO N=1,7  
        READ(1,*,IOSTAT=ISO) ISSPH(N),NPSPH(N),ISRSPH(N),ISPHXY(N)  
        WRITE(7,1002)NCARD  
        WRITE(7,*) ISSPH(N),NPSPH(N),ISRSPH(N),ISPHXY(N)  
      ENDDO  
      IF(ISO.GT.0) GOTO 100  
      ISSPH(8)=0  
      NPSPH(8)=1  
      DO N=1,7  
        IF(ISSPH(N).GE.1.AND.ISPHXY(N).EQ.3)ISSPH(8)=1  
      ENDDO  
      IF(ISSPH(8).EQ.1)THEN  
        DO N=1,7  
          NPSPH(8)=MAX(NPSPH(N),NPSPH(8))  
        ENDDO  
      ENDIF  
C  
C71A*  CONTROLS FOR HORIZONTAL PLANE SEDIMENT BED PROPERTIES
      NCARD='71A'  
      CALL SEEK('C71A')  
      READ(1,*,IOSTAT=ISO) ISBPH,ISBEXP,NPBPH,ISRBPH,ISBBDN,ISBLAY,  
     &    ISBPOR,ISBSED,ISBSND,ISBVDR,ISBARD  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ISBPH, ISBEXP,NPBPH,ISRBPH,ISBBDN,ISBLAY,  
     &    ISBPOR,ISBSED,ISBSND,ISBVDR,ISBARD  
      IF(ISO.GT.0) GOTO 100  
      IF(ISBEXP.GE.1) NPSPH(8)=MAX(NPSPH(8),NPBPH)  
      JSBPH=1  
      JSBPHA=1  
C  
C71B*  CONTROLS FOR FOOD CHAIN MODEL OUTPUT
      NCARD='71B'  
      CALL SEEK('C71B')  
      READ(1,*,IOSTAT=ISO) ISFDCH,NFDCHZ,HBFDCH,TFCAVG
      WRITE(7,1002)NCARD
      WRITE(7,*) ISFDCH,NFDCHZ,HBFDCH,TFCAVG
      IF(ISO.GT.0) GOTO 100
C
C72*  CONTROLS FOR HORIZONTAL PLANE SURFACE ELEVATION OR PRESSURE
	JSFDCH=1
      NCARD='72'  
      CALL SEEK('C72')  
      READ(1,*,IOSTAT=ISO) ISPPH,NPPPH,ISRPPH,IPPHXY  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ISPPH,NPPPH,ISRPPH,IPPHXY  
      IF(ISO.GT.0) GOTO 100  
C  
C73*  CONTROLS FOR HORIZONTAL PLANE VELOCITY PLOTTING 
      NCARD='73'  
      CALL SEEK('C73')  
      READ(1,*,IOSTAT=ISO) ISVPH,NPVPH,ISRVPH,IVPHXY  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ISVPH,NPVPH,ISRVPH,IVPHXY  
      IF(ISO.GT.0) GOTO 100  

      CALL CHECK_OUTPUT_FREQUENCY
#ifdef key_ncdf
      ! We need information on I & J location in cell to manage reconciliation
      ! Ensure these are present for velocity elevation and transported solutes
      DO N=1,7
        IF(ISSPH(N)==1.AND.ISPHXY(N)/=1)THEN
          IF(PARTID==MASTER_TASK)THEN
            PRINT*,'ISPHXY(',N,') on C71 SET TO 1'
          ENDIF
          ISPHXY(N)=1
        ENDIF
      ENDDO  
      IF(ISPPH==1.AND.IPPHXY/=1)THEN
        IF(PARTID==MASTER_TASK)THEN
          PRINT*,'IPPHXY on C72 SET TO 1'    
        ENDIF
        IPPHXY=1
      ENDIF 
      IF(ISVPH==1.AND.IVPHXY/=1)THEN
        IF(PARTID==MASTER_TASK)THEN
          PRINT*,'IVPHXY on C73 SET TO 1'    
        ENDIF
        IVPHXY=1
      ENDIF
#endif
C  
C74*  CONTROLS FOR VERTICAL PLANE SCALAR FIELD CONTOURING 
      NCARD='74'  
      CALL SEEK('C74')  
      READ(1,*,IOSTAT=ISO) ISECSPV,NPSPV(1),ISSPV(1),ISRSPV(1),  
     &    ISHPLTV(1)  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ISECSPV,NPSPV(1),ISSPV(1),ISRSPV(1),  
     &    ISHPLTV(1)  
      SHPLTV(1)=FLOAT(ISHPLTV(1))  
      SBPLTV(1)=1.0-SHPLTV(1)  
      DO N=2,7  
        READ(1,*,IOSTAT=ISO) IDUMMY,NPSPV(N),ISSPV(N),ISRSPV(N),  
     &      ISHPLTV(N)  
        WRITE(7,1002)NCARD  
        WRITE(7,*) IDUMMY,NPSPV(N),ISSPV(N),ISRSPV(N),  
     &      ISHPLTV(N)  
        SHPLTV(N)=FLOAT(ISHPLTV(N))  
        SBPLTV(N)=1.0-SHPLTV(N)  
      ENDDO  
      IF(ISO.GT.0) GOTO 100  
C  
      IF(ISECSPV.GT.0)THEN  
C75*  MORE CONTROLS FOR VERTICAL PLANE SCALAR FIELD CONTOURING
        NCARD='75'  
        CALL SEEK('C75')  
        DO IS=1,ISECSPV  
          READ(1,*,IOSTAT=ISO) DUM,NIJSPV(IS),CCTITLE(10+IS)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) DUM,NIJSPV(IS),CCTITLE(10+IS)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
C76*  I,J LOCATIONS DEFINING VERTICAL PLANE FOR CONTOURING 
        NCARD='76'  
        CALL SEEK('C76')  
        DO IS=1,ISECSPV  
          DO NPP=1,NIJSPV(IS)  
            READ(1,*,IOSTAT=ISO) DUM,ISPV(NPP,IS),JSPV(NPP,IS)  
            WRITE(7,1002)NCARD  
            WRITE(7,*) DUM,ISPV(NPP,IS),JSPV(NPP,IS)  
            IF(ISO.GT.0) GOTO 100  
          ENDDO  
        ENDDO  
      ENDIF  
C  
      NCARD='77'  
      CALL SEEK('C77')  
      READ(1,*,IOSTAT=ISO) ISECVPV,NPVPV,ISVPV,ISRVPV  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ISECVPV,NPVPV,ISVPV,ISRVPV  
      IF(ISO.GT.0) GOTO 100  
C  
      IF(ISECVPV.GT.0)THEN  
        NCARD='78'  
        CALL SEEK('C78')  
        DO IS=1,ISECVPV  
          READ(1,*,IOSTAT=ISO) DUM,NIJVPV(IS),ANGVPV(IS),CVTITLE(10+IS)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) DUM,NIJVPV(IS),ANGVPV(IS),CVTITLE(10+IS)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
        NCARD='79'  
        CALL SEEK('C79')  
        DO IS=1,ISECVPV  
          DO NPP=1,NIJVPV(IS)  
            READ(1,*,IOSTAT=ISO) DUM,IVPV(NPP,IS),JVPV(NPP,IS)  
            WRITE(7,1002)NCARD  
            WRITE(7,*) DUM,IVPV(NPP,IS),JVPV(NPP,IS)  
            IF(ISO.GT.0) GOTO 100  
          ENDDO  
        ENDDO  
      ENDIF  
C  
      NCARD='80'  
      CALL SEEK('C80')  
      READ(1,*,IOSTAT=ISO)IS3DO,ISR3DO,NP3DO,KPC,NWGG,I3DMIN,I3DMAX,  
     &    J3DMIN,J3DMAX,I3DRW,SELVMAX,BELVMIN  
      WRITE(7,1002)NCARD  
      WRITE(7,*)IS3DO,ISR3DO,NP3DO,KPC,NWGG,I3DMIN,I3DMAX,  
     &    J3DMIN,J3DMAX,I3DRW,SELVMAX,BELVMIN  
      IF(ISO.GT.0) GOTO 100  
      NCALL3D=0  
      NRCAL3D=0  
C  
      NCARD='81'  
      CALL SEEK('C81')  
      READ(1,*,IOSTAT=ISO)CDUM,IS3DUUU,JS3DUUU,UUU3DMA,UUU3DMI  
      WRITE(7,1002)NCARD  
      WRITE(7,*)CDUM,IS3DUUU,JS3DUUU,UUU3DMA,UUU3DMI  
      IF(ISO.GT.0) GOTO 100  
      READ(1,*,IOSTAT=ISO)CDUM,IS3DVVV,JS3DVVV,VVV3DMA,VVV3DMI  
      WRITE(7,1002)NCARD  
      WRITE(7,*)CDUM,IS3DVVV,JS3DVVV,VVV3DMA,VVV3DMI  
      IF(ISO.GT.0) GOTO 100  
      READ(1,*,IOSTAT=ISO)CDUM,IS3DWWW,JS3DWWW,WWW3DMA,WWW3DMI  
      WRITE(7,1002)NCARD  
      WRITE(7,*)CDUM,IS3DWWW,JS3DWWW,WWW3DMA,WWW3DMI  
      IF(ISO.GT.0) GOTO 100  
      READ(1,*,IOSTAT=ISO)CDUM,IS3DSAL,JS3DSAL,SAL3DMA,SAL3DMI  
      WRITE(7,1002)NCARD  
      WRITE(7,*)CDUM,IS3DSAL,JS3DSAL,SAL3DMA,SAL3DMI  
      IF(ISO.GT.0) GOTO 100  
      READ(1,*,IOSTAT=ISO)CDUM,IS3DTEM,JS3DTEM,TEM3DMA,TEM3DMI  
      WRITE(7,1002)NCARD  
      WRITE(7,*)CDUM,IS3DTEM,JS3DTEM,TEM3DMA,TEM3DMI  
      IF(ISO.GT.0) GOTO 100  
      READ(1,*,IOSTAT=ISO)CDUM,IS3DDYE,JS3DDYE,DYE3DMA,DYE3DMI  
      WRITE(7,1002)NCARD  
      WRITE(7,*)CDUM,IS3DDYE,JS3DDYE,DYE3DMA,DYE3DMI  
      IF(ISO.GT.0) GOTO 100  
      READ(1,*,IOSTAT=ISO)CDUM,IS3DSED,JS3DSED,SED3DMA,SED3DMI  
      WRITE(7,1002)NCARD  
      WRITE(7,*)CDUM,IS3DSED,JS3DSED,SED3DMA,SED3DMI  
      IF(ISO.GT.0) GOTO 100  
      READ(1,*,IOSTAT=ISO)CDUM,IS3DSND,JS3DSND,SND3DMA,SND3DMI  
      WRITE(7,1002)NCARD  
      WRITE(7,*)CDUM,IS3DSND,JS3DSND,SND3DMA,SND3DMI  
      IF(ISO.GT.0) GOTO 100  
      READ(1,*,IOSTAT=ISO)CDUM,IS3DTOX,JS3DTOX,TOX3DMA,TOX3DMI  
      WRITE(7,1002)NCARD  
      WRITE(7,*)CDUM,IS3DTOX,JS3DTOX,TOX3DMA,TOX3DMI  
      IF(ISO.GT.0) GOTO 100  
C  
      NCARD='82'  
      CALL SEEK('C82')  
      READ(1,*,IOSTAT=ISO) ISLSHA,MLLSHA,NTCLSHA,ISLSTR,ISHTA  
      WRITE(7,1002)NCARD  
      WRITE(7,*) ISLSHA,MLLSHA,NTCLSHA,ISLSTR,ISHTA  
      IF(ISO.GT.0) GOTO 100  
C  
      IF(MLLSHA.GT.0)THEN  
        NCARD='83'  
        CALL SEEK('C83')  
        DO M=1,MLLSHA  
          READ(1,*,IOSTAT=ISO) ILLSHA(M),JLLSHA(M),LSHAP(M),LSHAB(M),  
     &        LSHAUE(M),LSHAU(M),CLSL(M)  
          WRITE(7,1002)NCARD  
          WRITE(7,*) ILLSHA(M),JLLSHA(M),LSHAP(M),LSHAB(M),  
     &        LSHAUE(M),LSHAU(M),CLSL(M)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
      ENDIF  
C  
      NCARD='84'  
      CALL SEEK('C84')  
      READ(1,*,IOSTAT=ISO)ISTMSR,MLTMSR,NBTMSR,NSTMSR,NWTMSR,NTSSTSP,  
     &    TCTMSR  
      WRITE(7,1002)NCARD  
      WRITE(7,*)ISTMSR,MLTMSR,NBTMSR,NSTMSR,NWTMSR,NTSSTSP,  
     &    TCTMSR  
      IF(ISO.GT.0) GOTO 100
  
      JSTMSR=1  
      NCTMSR=1  
      JSHYDOUT=1  
      NCHYDOUT=1  
C  
      IF(NTSSTSP.GT.0)THEN  
        NCARD='85'  
        CALL SEEK('C85')  
        DO ITSSS=1,NTSSTSP  
          READ(1,*,IOSTAT=ISO)IDUM,MTSSTSP(ITSSS)  
          WRITE(7,1002)NCARD  
          WRITE(7,*)IDUM,MTSSTSP(ITSSS)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
C  
        NCARD='86'  
        CALL SEEK('C86')  
        DO ITSSS=1,NTSSTSP  
          DO MTSSS=1,MTSSTSP(ITSSS)  
            READ(1,*,IOSTAT=ISO)IDUM,IDUM,TSSTRT(MTSSS,ITSSS),  
     &          TSSTOP(MTSSS,ITSSS)  
            WRITE(7,1002)NCARD  
            WRITE(7,*)IDUM,IDUM,TSSTRT(MTSSS,ITSSS),  
     &          TSSTOP(MTSSS,ITSSS)  
            IF(ISO.GT.0) GOTO 100  
          ENDDO  
        ENDDO  
      ENDIF  
C  
      IF(MLTMSR.GT.0)THEN  
        NCARD='87'  
        CALL SEEK('C87')  
        DO M=1,MLTMSR  
      READ(1,*,IOSTAT=ISO)ILTMSR_GL(M),JLTMSR_GL(M),NTSSSS_GL(M),MTMSRP_GL(M),
     &       MTMSRC_GL(M),MTMSRA_GL(M),MTMSRUE_GL(M),MTMSRUT_GL(M),
     &       MTMSRU_GL(M),MTMSRQE_GL(M),MTMSRQ_GL(M),CLTMSR_GL(M)
      WRITE(7,1002)NCARD
      WRITE(7,*)ILTMSR_GL(M),JLTMSR_GL(M),NTSSSS_GL(M),MTMSRP_GL(M),
     &       MTMSRC_GL(M),MTMSRA_GL(M),MTMSRUE_GL(M),MTMSRUT_GL(M),
     &       MTMSRU_GL(M),MTMSRQE_GL(M),MTMSRQ_GL(M),CLTMSR_GL(M) 
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
      ENDIF  
C  
      NCARD='88'  
      CALL SEEK('C88')  
      READ(1,*,IOSTAT=ISO)ISVSFP,MDVSFP,MLVSFP,TMVSFP,TAVSFP  
      WRITE(7,1002)NCARD  
      WRITE(7,*)ISVSFP,MDVSFP,MLVSFP,TMVSFP,TAVSFP  
      IF(ISO.GT.0) GOTO 100  
      JSVSFP=1  
C  
      IF(MDVSFP.GT.0)THEN  
        NCARD='89'  
        CALL SEEK('C89')  
        DO M=1,MDVSFP  
          READ(1,*,IOSTAT=ISO)IDUM,DMVSFP(M)  
          WRITE(7,1002)NCARD  
          WRITE(7,*)IDUM,DMVSFP(M)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
      ENDIF  
C  
      IF(MLVSFP.GT.0)THEN  
        NCARD='90'  
        CALL SEEK('C90')  
        DO M=1,MLVSFP  
          READ(1,*,IOSTAT=ISO)IDUM,TIMVSFP(M),IVSFP(M),JVSFP(M)  
          WRITE(7,1002)NCARD  
          WRITE(7,*)IDUM,TIMVSFP(M),IVSFP(M),JVSFP(M)  
          IF(ISO.GT.0) GOTO 100  
        ENDDO  
      ENDIF  
C  
      NTS=NTC*NTSPTC  
      NBVSFP=NTC*NTSPTC  
      NSVSFP=0  
C  
      DO M=1,MLVSFP  
        TIMVSFP(M)=TMVSFP*(TIMVSFP(M)+TAVSFP)  
      ENDDO  
C  
      DT=TIDALP*FLOAT(NFLTMT)/FLOAT(NTSPTC)  
C  
      DO M=1,MLVSFP  
        NTMP=NINT( (TIMVSFP(M)-TCON*TBEGIN)/DT )  
        NTMP=MIN(NTMP,NTS)  
        NBVSFP=MIN(NBVSFP,NTMP)-1  
        NSVSFP=MAX(NSVSFP,NTMP)+1  
        NTVSFP(M)=NTMP  
      ENDDO  
C  
      DO M=1,MLVSFP  
        TIMVSFP(M)=(TIMVSFP(M)/TMVSFP)-TAVSFP  
      ENDDO  
      GOTO 2000  
C  
C **  WRITE INPUT ERROR MESSAGES AND TERMINATE RUN  
C  
  100 WRITE(6,1001)NCARD  
      WRITE(8,1001)NCARD  
      WRITE(7,1001)NCARD  
      STOP  
 2000 CONTINUE  
  
C **  NOW REWIND UNIT 1 & READ IN AS CHARACTER TO WRITE TO UNIT 7  
C  
      REWIND (1)  
   21 READ(1,22,END=24) TEXT  
      WRITE (7,23) TEXT  
      GOTO 21  
   24 CONTINUE  
      CLOSE(1)  
   22 FORMAT (A80)  
   23 FORMAT (1X,A80)  
C  
C **  READ CELL TYPES FROM FILES CELL.INP  
C  
      IF (PARTID ==0) PRINT *,'READING CELL.INP'
      OPEN(1,FILE='CELL.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES AND DETERMINE FILE FORMAT  
C  
      DO IS=1,4  
        READ(1,*)  
      ENDDO  
      READ(1,66)ADUMMY
      READ(ADUMMY,*)JCTMP
      
      CLOSE(1)  
      OPEN(1,FILE='CELL.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
      DO IS=1,4  
        READ(1,*)  
      ENDDO  
      IF(JCTMP.NE.JC_GLOBAL)THEN  
C  
C **    READ OLD FILE FORMAT  
C  
        JACROSS=JC_GLOBAL  
        IF(JC.GT.120)JACROSS=120  
        DO JT=1,JC,JACROSS  
          JF=JT  
          JLAST=JT+JACROSS-1  
          IF(JLAST.GT.JC) JLAST=JC  
          WRITE (7,8)JF,JLAST  
          DO I=1,IC  
            READ(1,6,IOSTAT=ISO) (IJCT(I,J),J=JF,JLAST)  
            IF(ISO.GT.0) GOTO 800  
            WRITE (7,16) (IJCT(I,J),J=JF,JLAST)  
          ENDDO  
          WRITE(7,15)  
        ENDDO  
      ELSE  
C  
        IF(IC_GLOBAL.GT.120)THEN  
          IACROSS=120  
          DO IT=1,IC_GLOBAL,IACROSS  
            IFIRST=IT  
            ILAST=IT+IACROSS-1  
            IF(ILAST.GT.IC_GLOBAL) ILAST=IC_GLOBAL  
            WRITE (7,88)IFIRST,ILAST  
            DO J=JC_GLOBAL,1,-1 
              READ(1,66,IOSTAT=ISO)ADUMMY,(IJCT_GLOBAL(I,J),I=IFIRST,ILAST)  
              IF(ISO.GT.0) GOTO 800  
              WRITE (7,166)ADUMMY,(IJCT_GLOBAL(I,J),I=IFIRST,ILAST)  
            ENDDO  
            WRITE(7,15)  
          ENDDO  
        ELSE  
          IFIRST=1  
          ILAST=IC_GLOBAL  
          WRITE (7,88)IFIRST,ILAST  
          DO J=JC_GLOBAL,1,-1  
            READ(1,66,IOSTAT=ISO)ADUMMY,(IJCT_GLOBAL(I,J),I=IFIRST,ILAST)  
            IF(ISO.GT.0) GOTO 800  
            WRITE (7,166)ADUMMY,(IJCT_GLOBAL(I,J),I=IFIRST,ILAST)  
          ENDDO  
          WRITE(7,15)  
        ENDIF  
      ENDIF  
      CLOSE(1)  
C  
    8 FORMAT ('   CELL TYPE ARRAY,J=',I5,2X,'TO J=',I5,//)  
C
C----------------------------------------------------------------------C
C
      OPEN(1,FILE='CELLLT.INP',STATUS='UNKNOWN')
      IF (PARTID ==0) PRINT *,'READING CELLLT.INP'

C
C **  SKIP OVER TITLE AND AND HEADER LINES AND DETERMINE FILE FORMAT
C 
      DO IS=1,4  
        READ(1,*)  
      ENDDO  
      READ(1,66)ADUMMY
      READ(ADUMMY,*)JCTMP

      CLOSE(1)
      OPEN(1,FILE='CELLLT.INP',STATUS='UNKNOWN')
C
C **  SKIP OVER TITLE AND AND HEADER LINES
C 
      DO IS=1,4
        READ(1,*)  
      ENDDO
C
      IF(JCTMP.NE.JC_GLOBAL)THEN
C **    READ OLD FILE FORMAT
        JACROSS=JC
        IF(JC.GT.120)JACROSS=120
        DO JT=1,JC,JACROSS
          JF=JT
          JLAST=JT+JACROSS-1
          IF(JLAST.GT.JC) JLAST=JC
          WRITE (7,8)JF,JLAST 
          DO I=1,IC
            READ(1,6,IOSTAT=ISO) (IJCTLT(I,J),J=JF,JLAST)
            IF(ISO.GT.0) GOTO 800
            WRITE (7,16) (IJCTLT(I,J),J=JF,JLAST)
          ENDDO
          WRITE(7,15)
        ENDDO

      ELSE
C **    READ NEW FILE FORMAT
        IF(IC_GLOBAL.GT.120)THEN
          IACROSS=120
          DO IT=1,IC_GLOBAL,IACROSS
            IFIRST=IT
            ILAST=IT+IACROSS-1
            IF(ILAST.GT.IC_GLOBAL) ILAST=IC_GLOBAL
            WRITE (7,88)IFIRST,ILAST 
            DO J=JC_GLOBAL,1,-1
              READ(1,66,IOSTAT=ISO)ADUMMY,(IJCTLT_GLOBAL(I,J),I=IFIRST,ILAST)
              IF(ISO.GT.0) GOTO 800
              WRITE (7,166)ADUMMY,(IJCTLT_GLOBAL(I,J),I=IFIRST,ILAST)
            ENDDO
            WRITE(7,15)
          ENDDO
        ELSE
          IFIRST=1
          ILAST=IC_GLOBAL
          WRITE (7,88)IFIRST,ILAST 
          DO J=JC_GLOBAL,1,-1
            READ(1,66,IOSTAT=ISO)ADUMMY,(IJCTLT_GLOBAL(I,J),I=IFIRST,ILAST)
            IF(ISO.GT.0) GOTO 800
            WRITE (7,166)ADUMMY,(IJCTLT_GLOBAL(I,J),I=IFIRST,ILAST)
          ENDDO
          WRITE(7,15)
        ENDIF
      ENDIF
C
      CLOSE(1)
C
   88 FORMAT ('   CELLLT TYPE ARRAY,I=',I5,2X,'TO I=',I5,//)  
C  
C **  IF ISPGNS GE 1, READ IN NORTH-SOUTH BOUNDARY CELLS FROM  
C **  FILE MAPPGNS.INP TO SPECIFY A PERIODIC DOMAIN IN THE NORTH-SOUTH  
C **  DIRECTION  
C  
      IF(ISPGNS.GE.1)THEN  
        PRINT *,'READING MAPPGNS.INP'
        OPEN(1,FILE='MAPPGNS.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,8  
          READ(1,*)  
        ENDDO  
        READ(1,*,IOSTAT=ISO) NPNSBP  
        IF(ISO.GT.0) GOTO 950  
        DO NPNS=1,NPNSBP  
          READ(1,*,IOSTAT=ISO) ISPNS(NPNS),JSPNS(NPNS),  
     &        INPNS(NPNS),JNPNS(NPNS)  
          IF(ISO.GT.0) GOTO 950  
        ENDDO  
        CLOSE(1)  
      ENDIF  
C  
C **  GENERATE CELL MAPPINGS  
C  
      CALL CELLMAP  


      IF(NASER.LE.1)THEN  
        DO L=2,LA  
          ATMWHT(L,1)=1.  
        ENDDO  
      ENDIF  
      IF(NASER.GT.1)THEN  
        PRINT *,'READING ATMMAP.INP'
        OPEN(121,FILE='ATMMAP.INP',STATUS='UNKNOWN')  
        DO IS=1,4  
          READ(121,*)  
        ENDDO  
        DO LGLOB=2,LC_GLOBAL-1 
          READ(121,*)LD,ID,JD,(ATMWHT_TEMP(N),N=1,NASER) 
          I = XLOC(ID)
          J = YLOC(JD)
          IF (I > 0 .AND. I <=IC) THEN
            IF (J > 0 .AND.J <=JC) THEN
              L = LIJ(I,J)
              ATMWHT(L,1:KC) = ATMWHT_TEMP(1:KC)
            END IF
          END IF
        ENDDO  
        CLOSE(121)  
      ENDIF
!      DO L = 2,LA
!         IF (PARTID == 1) WRITE(*,*)XPAR(IL(L)), YPAR(JL(L)), SUM((ATMWHT(L,:)))
!      ENDDO 
      IF(NWSER<=1)THEN  
        DO L=2,LA  
          WNDWHT(L,1)=1.  
        ENDDO  
      ENDIF  
      IF(NWSER.GT.1)THEN  
        PRINT *,'READING WNDMAP.INP'
        OPEN(122,FILE='WNDMAP.INP',STATUS='UNKNOWN')  
        DO IS=1,4  
          READ(122,*)  
        ENDDO 
         DO LGLOB=2,LC_GLOBAL-1 
          READ(122,*)LD,ID,JD,(WNDWHT_TEMP(N),N=1,NWSER) 
          I = XLOC(ID)
          J = YLOC(JD)
          IF (I > 0 .AND. I <=IC) THEN
             IF (J > 0 .AND.J <=JC) THEN
               L = LIJ(I,J)
               WNDWHT(L,1:KC) = WNDWHT_TEMP(1:KC)
            END IF
          END IF
        ENDDO  
      CLOSE(122)  
      ENDIF  
 

C     DHC
      IF(ISWAVE.EQ.3.OR.ISPD.GE.2)THEN
        ! *** COMPUTE CELL AREAS AND CENTROIDS, REQUIRES CORNERS.INP FILE
        CALL AREA_CENTRD
      ENDIF
      
C  
C **  READ IN CELL CENTER DEPTHS AND BOTTOM BED ELEVATION FOR CARTESIAN  
C **  OR MIXED CATERSIAN CURVILINEAR GRID FORM FILE DEPTH.INP.  IF THE  
C **  NUMBER OF VARIALBLE CELLS EQUALS NUMBER OF WATER CELLS OR  
C **  ISCLO=1 DEPTHS ARE READ FROM THE FILE DXDY.INP  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C       READ(1,7,IOSTAT=ISO) (RTMP1IJ(I,J),J=1,JC)  
C      ELSE  
C       READ(1,7,IOSTAT=ISO) (RTMP1IJ(I,J),J=JF,JLAST)  
C **  SKIP OVER BLANK, TITLE AND AND HEADER LINES  
C       READ(1,7,IOSTAT=ISO) (RTMP1IJ(I,J),J=1,JC)  
C      ELSE  
C       READ(1,7,IOSTAT=ISO) (RTMP1IJ(I,J),J=JF,JLAST)  
C **  SKIP OVER BLANK, TITLE AND AND HEADER LINES  
C       READ(1,7,IOSTAT=ISO) (RTMP1IJ(I,J),J=1,JC)  
C      ELSE  
C       READ(1,7,IOSTAT=ISO) (RTMP1IJ(I,J),J=JF,JLAST)  
C  
   15 FORMAT (/)  
    6 FORMAT (120I1)  
C   66 FORMAT (I3,2X,120I1)     PMC
   66 FORMAT (A5,120I1)  
    9 FORMAT (/,' DEPTH ARRAY:',//)  
   16 FORMAT (1X,120I1)  
C  166 FORMAT (1X,I3,2X,120I1)  PMC
  166 FORMAT (1X,A5,120I1)  
    7 FORMAT (30F4.1)  
   17 FORMAT(1X,30F4.1)  
C  
C **  READ CURVILINEAR-ORTHOGONAL OR VARIABLE CELL DATA FROM FILE  
C **  DXDY.INP  
C **  INITIALIZE CELL DIMENSIONS TO CONSTANT CARTESIAN OR DUMMY VALUES  
C  
      DO L=1,LC  
        DXP(L)=DX*DXYCVT  
        DYP(L)=DY*DXYCVT  
        ZBR(L)=ZBRADJ  
      ENDDO  
C  
C **  READ IN DX, DY, DEPTH AND BOTTOM ELEVATION AT CELL CENTERS OF  
C **  VARIABLE CELLS  
C  
      LMHK=.FALSE.
      IF(LVC.GT.0)THEN 
        IF (PARTID ==0) PRINT *,'READING DXDY.INP'
        OPEN(1,FILE='DXDY.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C 
        DXP_GLOB(0:GNX,0:GNY)=1.0; DYP_GLOB(0:GNX,0:GNY)=1.0 !Initialize vectors
        HMP_GLOB(0:GNX,0:GNY)=HMIN !Initialize matrix
        DO IS=1,4  
          READ(1,*)  
        ENDDO  
        DO LT=1,LVC 
          IF(ISVEG.EQ.0)THEN  
            READ(1,*,IOSTAT=ISO)IGL,JGL,DXIJ,DYIJ,HIJ,BELVIJ,ZBRIJ
            IF(ISO.GT.0) GOTO 830
          ELSE
            READ(1,*,IOSTAT=ISO)IGL,JGL,DXIJ,DYIJ,HIJ,BELVIJ,ZBRIJ,MVEGIJT !SCJ adding MHK devices when MVEGIJT>90  
            IF(ISO.GT.0) GOTO 830
          ENDIF
          I = XLOC(IGL)
          J = YLOC(JGL)
          IF (I >0 .AND. I <= IC) THEN
            IF (J > 0 .AND. J <= JC )THEN  
              L=LIJ(I,J) 
              IF (L > 0) THEN
                DXP(L)=DXYCVT*DXIJ  
                DYP(L)=DXYCVT*DYIJ  
                HMP(L)=HADADJ + HCVRT*HIJ  
                HMP(L)=MAX(HMP(L),HMIN)  
                DXP_GLOB(IGL,JGL) = DXYCVT*DXIJ
                DYP_GLOB(IGL,JGL) = DXYCVT*DYIJ
                HMP_GLOB(IGL,JGL) = HMP(L)
                BELV(L)=BELADJ + BELCVRT*BELVIJ  
                BELV1(L)=BELADJ + BELCVRT*BELVIJ  
                ZBR(L)=ZBRADJ + ZBRCVRT*ZBRIJ
                IF(ISVEG==1)THEN
                  MVEGL(L)=MVEGIJT
		        IF(MVEGIJT>90)THEN !MHK IF STATEMENT
                    LMHK=.TRUE.
                    ITURB=ITURB+1
		          IJLTURB(ITURB,1)=I
		          IJLTURB(ITURB,2)=J
		          IJLTURB(ITURB,3)=LIJ(I,J)
                  ENDIF
                ENDIF
              ENDIF
            ENDIF
          ENDIF
        ENDDO 
      ENDIF
      DO I = 1,IC_GLOBAL
        DO J = 1,JC_GLOBAL
        ILOC = XLOC(I)
        JLOC = YLOC(J)
        IF (ILOC .GT. 0 .AND. ILOC .LE. IC) THEN
          IF (JLOC .GT. 0 .AND. JLOC .LE. JC) THEN
            L = LIJ(ILOC,JLOC)
            IF ( L > 0 ) THEN
              LS=LSC(L)  
              DXU(L)=0.5*(DXP_GLOB(I,J)+DXP_GLOB(I-1,J))  
              DYU(L)=0.5*(DYP_GLOB(I,J)+DYP_GLOB(I-1,J))  
              DXV(L)=0.5*(DXP_GLOB(I,J)+DXP_GLOB(I,J-1))  
              DYV(L)=0.5*(DYP_GLOB(I,J)+DYP_GLOB(I,J-1))  
              HMU(L)=0.5*(DXP(L)*DYP(L)*HMP(L)+DXP_GLOB(I-1,J)*DYP_GLOB(I-1,J)*  
     &            HMP_GLOB(I-1,J))/(DXU(L)*DYU(L))  
              HMV(L)=0.5*(DXP(L)*DYP(L)*HMP(L)+DXP_GLOB(I,J-1 )*DYP_GLOB(I,J-1 )*HMP_GLOB(I,J-1 ))  
     &            /(DXV(L)*DYV(L))  
            END IF
          END IF
        END IF
        ENDDO  
      END DO

C#######################################################################
C     SEI Change, 1204, CAJ SEDZLJ
C     Change to implement Sedflume Data
C     Located here so spatial data can be incorporated after
C     the grid is set up
C     READ IN Sedflume Initial Data
C     FLAG set for Sedflume or not in SEDIC
C
      NSEDFLUME=0
      IF (IWRSP(1).EQ.98) THEN
	   NSEDFLUME=1
	   CALL SEDIC
	ELSEIF (IWRSP(1).EQ.99) THEN
	   NSEDFLUME=2
	   CALL SEDIC
	ENDIF
!***************************************************************
!PT check the SED distribution in layers for particle class 1. 7/24/08
      IF (NSEDFLUME.GT.0)THEN
        OPEN(UNIT=149,FILE='sedcheck.dat')
        WRITE(149,*)'size class', 1
        WRITE(149,*)'time    ','waterheight   ','  expval  ','DWS',
     &'1stlayer   ','2ndlayer   ','3rdlayer   ','4thlayer  ','5thlayer',
     &  'waterwieght1','waterwieght2','waterwieght3','waterwieght4',
     &  'waterwieght5','shear'
!***************************************************************
        inquire(file='ensight.inp',exist=status)
        if(status) then
        ! PT create ensight output files.
           OPEN(UNIT=117,FILE='ensight.inp',STATUS='old')
           ISSKIP = 0  
           READ(117,'(A1)') CCMRM  
           BACKSPACE(1)  
           IF(CCMRM .EQ. '#') ISSKIP = 1   
           !IF(ISSKIP .GT. 0) CALL SKIPCOMM(117,CCMRM)
           IF(ISSKIP .GT. 0) CALL SKIPCOMM(117,'#')
           READ(117,*,IOSTAT=ISO)ENSIGHT1, ENSIGHT2, ENSIGHT3, ENSIGHT4,
     & ENSIGHT5,  ENSIGHT6 , ENSIGHT7,  ENSIGHT8,  ENSIGHT9,  ENSIGHT10,
     & ENSIGHT11, ENSIGHT12, ENSIGHT13, ENSIGHT14, ENSIGHT15, ENSIGHT16,
     & ENSIGHT17, ENSIGHT18, ENSIGHT19, ENSIGHT20, ENSIGHT21, ENSIGHT22,
     & ENSIGHT23, ENSIGHT24, ENSIGHT25, ENSIGHT26, ENSIGHT27, ENSIGHT28,
     & ENSIGHT29, ENSIGHT30, ENSIGHT31
           CLOSE(117)
        
           OPEN(UNIT=113,FILE='ensight.case',FORM='FORMATTED')
           OPEN(UNIT=114,FILE='ensight.geo',FORM='FORMATTED')
           !******VARIABLES FILES*******
           IF(ENSIGHT1.GT.0) THEN
              OPEN(UNIT=118,FILE='ensight.U',FORM='FORMATTED')
           ENDIF
           IF(ENSIGHT2.GT.0) THEN
              OPEN(UNIT=119,FILE='ensight.V',FORM='FORMATTED')
           ENDIF
           IF( ISTRAN(6) .GT. 0) THEN
              IF(ENSIGHT3.GT.0) THEN
                 OPEN(UNIT=120,FILE='ensight.TAU',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT4.GT.0) THEN
                 OPEN(UNIT=121,FILE='ensight.D50',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT5.GT.0) THEN
                 OPEN(UNIT=122,FILE='ensight.CBL',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT6.GT.0) THEN
                 OPEN(UNIT=123,FILE='ensight.SED',FORM='FORMATTED')
              ENDIF
           ENDIF
           IF(ISTRAN(8) .GT. 0) THEN
              IF(ENSIGHT7.GT.0) THEN
                 OPEN(UNIT=124,FILE='ensight.CHC',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT8.GT.0) THEN
                 OPEN(UNIT=125,FILE='ensight.CHD',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT9.GT.0) THEN
                 OPEN(UNIT=126,FILE='ensight.CHG',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT10.GT.0) THEN
                 OPEN(UNIT=127,FILE='ensight.ROC',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT11.GT.0) THEN
                 OPEN(UNIT=128,FILE='ensight.LOC',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT12.GT.0) THEN
                 OPEN(UNIT=129,FILE='ensight.DOC',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT13.GT.0) THEN
                 OPEN(UNIT=130,FILE='ensight.ROP',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT14.GT.0) THEN
                 OPEN(UNIT=131,FILE='ensight.LOP',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT15.GT.0) THEN
                 OPEN(UNIT=132,FILE='ensight.DOP',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT16.GT.0) THEN
                 OPEN(UNIT=133,FILE='ensight.P40',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT17.GT.0) THEN
                 OPEN(UNIT=134,FILE='ensight.RON',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT18.GT.0) THEN
                 OPEN(UNIT=135,FILE='ensight.LON',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT19.GT.0) THEN
                 OPEN(UNIT=136,FILE='ensight.DON',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT20.GT.0) THEN
                 OPEN(UNIT=137,FILE='ensight.NHX',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT21.GT.0) THEN
                 OPEN(UNIT=138,FILE='ensight.NOX',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT22.GT.0) THEN
                 OPEN(UNIT=139,FILE='ensight.SUU',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT23.GT.0) THEN
                 OPEN(UNIT=140,FILE='ensight.SAA',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT24.GT.0) THEN
                 OPEN(UNIT=141,FILE='ensight.COD',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT25.GT.0) THEN
                 OPEN(UNIT=142,FILE='ensight.DOX',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT26.GT.0) THEN
                 OPEN(UNIT=143,FILE='ensight.TM',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT27.GT.0) THEN
                 OPEN(UNIT=144,FILE='ensight.SCB',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT28.GT.0) THEN
                 OPEN(UNIT=145,FILE='ensight.CO2',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT29.GT.0) THEN
                 OPEN(UNIT=146,FILE='ensight.MAC',FORM='FORMATTED')
              ENDIF
              IF(ENSIGHT30.GT.0) THEN
                 OPEN(UNIT=147,FILE='ensight.HEAT',FORM='FORMATTED')
              ENDIF
           ENDIF
           IF( ISTRAN(2) .GT. 0) THEN
              IF(ENSIGHT31 .GT. 0) THEN
                 OPEN(UNIT=148,FILE='ensight.TEMP',FORM='FORMATTED')
              ENDIF
           ENDIF 
        ENDIF  ! *** STATUS
      ENDIF    ! *** NSEDFLUME 
C  
C **  OPEN FILE MODDXDY.INP TO MODIFY INPUT VALUES OF DX AND DY  
C  
      IF(IMDXDY.GT.0)THEN  
        PRINT *,'READING MODDXDY.INP'
        OPEN(1,FILE='MODDXDY.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND HEADER LINES  
C  
        DO IS=1,4  
          READ(1,*)  
        ENDDO  
        READ(1,*) NMDXDY  
        IF(NMDXDY.GE.1)THEN  
          DO NMD=1,NMDXDY  
            READ(1,*)ITMP,JTMP,RMDX,RMDY  
            LTMP=LIJ(ITMP,JTMP)  
            DXP(LTMP)=RMDX*DXP(LTMP)  
            DYP(LTMP)=RMDY*DYP(LTMP)  
          ENDDO  
        ENDIF  
        CLOSE(1)  
      ENDIF  
C  
C **  OPEN FILE MODCHAN.INP TO INSERT SUBGRID CHANNELS INTO  
C **  HOST CELLS  
C  
      MDCHH=0  
      IF(ISCHAN.GT.0)THEN  
        PRINT *,'READING MODCHAN.INP'
        OPEN(1,FILE='MODCHAN.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,8  
          READ(1,*)  
        ENDDO  
        IF(ISCHAN.EQ.1)THEN  
          READ(1,*) MDCHH,MDCHHD,MDCHHD2  
          READ(1,*) MDCHITM,MDCHHQ,QCHERR  
          IF(MDCHH.GE.1)THEN  
            DO NMD=1,MDCHH  
              READ(1,*)MDCHTYP(NMD),IMDCHH(NMD),JMDCHH(NMD),  
     &            IMDCHU(NMD),JMDCHU(NMD),  
     &            IMDCHV(NMD),JMDCHV(NMD)  
              QCHANU(NMD)=0.  
              QCHANUN(NMD)=0.  
              QCHANV(NMD)=0.  
              QCHANVN(NMD)=0.  
            ENDDO  
          ENDIF  
        ENDIF  
        IF(ISCHAN.EQ.2)THEN  
          READ(1,*) MDCHH,MDCHHD,MDCHHD2  
          READ(1,*) MDCHITM,MDCHHQ,QCHERR  
          IF(MDCHH.GE.1)THEN  
            DO NMD=1,MDCHH  
              READ(1,*)MDCHTYP(NMD),IMDCHH(NMD),JMDCHH(NMD),  
     &            IMDCHU(NMD),JMDCHU(NMD),  
     &            IMDCHV(NMD),JMDCHV(NMD),  
     &            CHANLEN(NMD),PMDCH(NMD)  
              QCHANU(NMD)=0.  
              QCHANUN(NMD)=0.  
              QCHANV(NMD)=0.  
              QCHANVN(NMD)=0.  
            ENDDO  
          ENDIF  
        ENDIF  
        CLOSE(1)  
        IF(MDCHH.GE.1)THEN  
          DO NMD=1,MDCHH  
            LMDCHH(NMD)=LIJ(IMDCHH(NMD),JMDCHH(NMD))  
            IF(IMDCHU(NMD).EQ.1.AND.JMDCHU(NMD).EQ.1)THEN  
              LMDCHU(NMD)=1  
            ELSE  
              LMDCHU(NMD)=LIJ(IMDCHU(NMD),JMDCHU(NMD))  
            ENDIF  
            IF(IMDCHV(NMD).EQ.1.AND.JMDCHV(NMD).EQ.1)THEN  
              LMDCHV(NMD)=1  
            ELSE  
              LMDCHV(NMD)=LIJ(IMDCHV(NMD),JMDCHV(NMD))  
            ENDIF  
          ENDDO  
        ENDIF  
      ENDIF  
C  
C **  OPEN FILE CHANSEC.INP FOR 1-D CHANNEL CROSS SECTION DATA  
C
C *** REMOVED 2004-09-19  PMC  
C  
C **  OPEN FILE GWATER.INP TO SPECIFY GROUNDWATER INTERACTION  
C **  BY INFILTRATION AND EVAPOTRANSPIRATION  
C  
      ISGWIE=0  
      IF(ISGWIT.EQ.1)THEN  
        PRINT *,'READING GWATER.INP'
        OPEN(1,FILE='GWATER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,6  
          READ(1,*)  
        ENDDO  
        READ(1,*) ISGWIE  
        IF(ISGWIE.GE.1)THEN  
          READ(1,*) DAGWZ,RNPOR,RIFTRM  
        ELSE  
          DAGWZ=0.0  
          RNPOR=1.E-12  
          RIFTRM=0.0  
        ENDIF  
        CLOSE(1)  
      ENDIF  
  339 FORMAT(2I5,6F14.5)  
C  
C **  OPEN FILE FBODY.INP TO READ IN SPATIALLY VARYING BODY FORCES  
C  
      IF(ISBODYF.GE.1)THEN  
        PRINT *,'READING FBODY.INP'
        OPEN(1,FILE='FBODY.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,7  
          READ(1,*)  
        ENDDO  
        READ(1,*)NBDF,CVTFACX,CVTFACY
        DO LL=1,NBDF  
          READ(1,*)ITMP,JTMP,FBODY1,FBODY2  
          L=LIJ(ITMP,JTMP)  
          FBODYFX(L,1:KC)=CVTFACX*FBODY1  
          FBODYFY(L,1:KC)=CVTFACY*FBODY2  
        END DO  
        FBODYFX(1,1:KC)=0.0  
        FBODYFY(1,1:KC)=0.0  
        FBODYFX(LC,1:KC)=0.0  
        FBODYFY(LC,1:KC)=0.0  
        CLOSE(1)  
      ELSE  
        DO L=1,LC  
          FBODYFX(L,1:KC)=0.0  
          FBODYFY(L,1:KC)=0.0  
        END DO  
      ENDIF  
C  
C **  OPEN FILE SEDBLBC.INP TO READ IN SEDIMENT BEDLOAD OUTFLOW  
C **  OR RECIRCULATION BOUNDARY CONDITIONS  
C  
      NSBDLDBC=0  
      IF(ISBDLDBC.GE.1)THEN  
        PRINT *,'READING SEDBLBC.INP'
        OPEN(1,FILE='SEDBLBC.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,12  
          READ(1,*)  
        ENDDO  
        READ(1,*)NSBDLDBC  
        DO N=1,NSBDLDBC  
          READ(1,*)ITMPU,JTMPU,ITMPD,JTMPD,ISDBLDIR(N)  
          LSBLBCU(N)=LIJ(ITMPU,JTMPU)  
          IF(ITMPD.GT.0.AND.JTMPD.GT.0) THEN  
            LSBLBCD(N)=LIJ(ITMPD,JTMPD)  
          ELSE  
            LSBLBCD(N)=0  
          ENDIF  
        ENDDO  
        CLOSE(1)  
      ENDIF  
C  
C **  OPEN FILE GWMAP.INP TO SPECIFY GROUNDWATER INTERACTION BY  
C **  AMBIENT GROUNDWATER FLOW  
      ! *** DSLLC Begin
      IF(ISGWIT.GT.1)THEN  
        IF(ISGWIT.EQ.2)THEN
          ISGWIE=2
        ELSE
          ISGWIE=0
        ENDIF

        IF(ISGWIT.EQ.3)THEN
          OPEN(1,FILE='GWSEEP.INP',STATUS='UNKNOWN')

          DO IS=1,4
            READ(1,2) TEXT
          ENDDO

          READ(1,*)NSEEPCLASSES
          DO M=1,NSEEPCLASSES  
            READ(1,*)IDUM,SEEPRATE(M)
          ENDDO  
          CLOSE(1)
        ENDIF
        PRINT *,'READING GWMAP.INP'
        OPEN(1,FILE='GWMAP.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,10  
          READ(1,*)  
        ENDDO  
        DO LL=2,LA  
          READ(1,*)ITMP,JTMP,IVALUE,RVALUE  
          L=LIJ(ITMP,JTMP)  
          IF(ISGWIT.EQ.3)THEN
            IF(IVALUE.GT.NSEEPCLASSES)THEN
              WRITE(6,*)'BAD SEEPAGE CLASS AT I,J=',ITMP,JTMP
              STOP
            ENDIF
            RIFTR(L)=SEEPRATE(IVALUE)*RVALUE
          ELSE
            NGWSL(L)=IVALUE  
            GWFAC(L)=RVALUE  
          ENDIF
        END DO  
        CLOSE(1)  
      ENDIF  
C
C **  READ IN SPATIALLY VARYING SEDIMENT ROUGHNESS HEIGHT FOR 
C **  DETERMINING GRAIN STRESS
C
      IF(ISTRAN(6).GE.1.OR.ISTRAN(7).GE.1)THEN
        IF(ISBEDSTR.EQ.3)THEN
          PRINT *,'READING SEDROUGH.INP'
          OPEN(1,FILE='SEDROUGH.INP')
          DO IS=1,2
            READ(1,*)
          ENDDO
          DO L=2,LC-1
            READ(1,*) LDUM,IDUM,JDUM,ZBRSED(L)
            IF(ZBRSED(L).LE.0.0)THEN
              STOP ' BAD SEDIMENT ROUGHNESS IN SEDROUGH.INP'
            ENDIF
          ENDDO
          CLOSE(1)
        ENDIF
      ENDIF
C
C  
C **  OPEN FILE DOCW.INP TO SPECIFY SPATIAL VARYING, TIME CONSTANT  
C **  DISSOLVED ORGANIC CARBON IN WATER COLUMN  
C  
      IVAL=0  
      DO NT=1,NTOX  
        IF(ISTOC(NT).EQ.1.OR.ISTOC(NT).EQ.2)IVAL=1  
      ENDDO  
      IF(IVAL.EQ.1)THEN  
        IF(ISTDOCW.EQ.1)THEN  
          PRINT *,'READING DOCW.INP'
          OPEN(1,FILE='DOCW.INP',STATUS='UNKNOWN')  
          DO IS=1,8  
            READ(1,*)  
          ENDDO  
          READ(1,*)ISALTYP  
          IF(ISALTYP.EQ.0)THEN  
            DO L=2,LC-1  
              READ(1,*,IOSTAT=ISO) (STDOCW(L,K),K=1,KC)  
              IF(ISO.GT.0) GOTO 854  
            ENDDO  
          ELSE  
            DO L=2,LC-1  
              READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(STDOCW(L,K),K=1,KC)  
              IF(ISO.GT.0) GOTO 854  
            ENDDO  
          ENDIF  
        ENDIF  
      ENDIF  
C  
C **  OPEN FILE POCW.INP TO SPECIFY SPATIAL VARYING, TIME CONSTANT  
C **  PARTICULATE ORGANIC CARBON IN WATER COLUMN  
C  
      IVAL=0  
      DO NT=1,NTOX  
        IF(ISTOC(NT).EQ.1)IVAL=1  
      ENDDO  
      IF(IVAL.EQ.1)THEN  
        IF(ISTPOCW.EQ.1)THEN  
          PRINT *,'READING POCW.INP'
          OPEN(1,FILE='POCW.INP',STATUS='UNKNOWN')  
          DO IS=1,8  
            READ(1,*)  
          ENDDO  
          READ(1,*)ISALTYP  
          IF(ISALTYP.EQ.0)THEN  
            DO L=2,LC-1  
              READ(1,*,IOSTAT=ISO) (STPOCW(L,K),K=1,KC)  
              IF(ISO.GT.0) GOTO 854  
            ENDDO  
          ELSE  
            DO L=2,LC-1  
              READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(STPOCW(L,K),K=1,KC)  
              IF(ISO.GT.0) GOTO 854  
            ENDDO  
          ENDIF  
        ENDIF  
      ENDIF  
C  
C **  OPEN FILE FPOCW.INP TO SPECIFY SPATIAL VARYING, TIME CONSTANT  
C **  PARTICULATE ORGANIC CARBON FRACTION FOR EACH SEDIMENT CLASS  
C **  IN WATER COLUMN  
C  
      IVAL=0  
      DO NT=1,NTOX  
        IF(ISTOC(NT).EQ.2.OR.ISTOC(NT).EQ.3)IVAL=1  
      ENDDO  
      IF(IVAL.EQ.1)THEN  
        IF(ISTPOCW.EQ.3)THEN  
          PRINT *,'READING FPOCW.INP'
          OPEN(1,FILE='FPOCW.INP',STATUS='UNKNOWN')  
          DO NS=1,NSED+NSND  
            DO IS=1,8  
              READ(1,*)  
            ENDDO  
            READ(1,*)ISALTYP  
            IF(ISALTYP.EQ.0)THEN  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO) (STFPOCW(L,K,NS),K=1,KC)  
                IF(ISO.GT.0) GOTO 854  
              ENDDO  
            ELSE  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &              (STFPOCW(L,K,NS),K=1,KC)  
                IF(ISO.GT.0) GOTO 854  
              ENDDO  
            ENDIF  
          ENDDO  
        ENDIF  
      ENDIF  
C  
C **  OPEN FILE DOCB.INP TO SPECIFY SPATIAL VARYING, TIME CONSTANT  
C **  DISSOLVED ORGANIC CARBON IN SEDIMENT BED  
C  
      IVAL=0  
      DO NT=1,NTOX  
        IF(ISTOC(NT).EQ.1.OR.ISTOC(NT).EQ.2)IVAL=1  
      ENDDO  
      IF(IVAL.EQ.1)THEN  
        IF(ISTDOCB.EQ.1)THEN  
          PRINT *,'READING DOCB.INP'
          OPEN(1,FILE='DOCB.INP',STATUS='UNKNOWN')  
          DO IS=1,8  
            READ(1,*)  
          ENDDO  
          READ(1,*)ISALTYP,IREAD,KBINPUT  
          DO K=1,KB  
            DO L=2,LC-1  
              STDOCB(L,K)=0.0  
            ENDDO  
          ENDDO  
          IF(IREAD.EQ.0)THEN  
            IF(ISALTYP.EQ.0)THEN  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO) STDOCB(L,1)  
                IF(ISO.GT.0) GOTO 856  
                DO K=2,KB  
                  STDOCB(L,K)=STDOCB(L,1)  
                ENDDO  
              ENDDO  
            ELSE  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,STDOCB(L,1)  
                IF(ISO.GT.0) GOTO 856  
                DO K=2,KB  
                  STDOCB(L,K)=STDOCB(L,1)  
                ENDDO  
              ENDDO  
            ENDIF  
          ENDIF  
          IF(IREAD.EQ.1)THEN  
            IF(ISALTYP.EQ.0)THEN  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO) (STDOCB(L,K),K=1,KB)  
                IF(ISO.GT.0) GOTO 856  
              ENDDO  
            ELSE  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &              (STDOCB(L,K),K=1,KB)  
                IF(ISO.GT.0) GOTO 856  
              ENDDO  
            ENDIF  
          ENDIF  
          IF(IREAD.EQ.2)THEN  
            IF(ISALTYP.EQ.0)THEN  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO) (STDOCB(L,K),K=1,KBINPUT)  
                IF(ISO.GT.0) GOTO 856  
                DO K=KBINPUT,KB  
                  STDOCB(L,K)=STDOCB(L,KBINPUT)  
                ENDDO  
              ENDDO  
            ELSE  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &              (STDOCB(L,K),K=1,KBINPUT)  
                IF(ISO.GT.0) GOTO 856  
                DO K=KBINPUT,KB  
                  STDOCB(L,K)=STDOCB(L,KBINPUT)  
                ENDDO  
              ENDDO  
            ENDIF  
          ENDIF  
          CLOSE(1)  
        ENDIF  
      ENDIF  
C  
C **  OPEN FILE POCB.INP TO READ SPATIALY VARYING, TIME CONSTANT  
C **  PARTICULATE ORGANIC CARBON IN BED  
C  
      IVAL=0  
      DO NT=1,NTOX  
        IF(ISTOC(NT).EQ.1)IVAL=1  
      ENDDO  
      IF(IVAL.EQ.1)THEN  
        IF(ISTPOCB.EQ.1)THEN  
          PRINT *,'READING POCB.INP'
          OPEN(1,FILE='POCB.INP',STATUS='UNKNOWN')  
          DO IS=1,8  
            READ(1,*)  
          ENDDO  
          READ(1,*)ISALTYP,IREAD,KBINPUT  
          DO K=1,KB  
            DO L=2,LC-1  
              STPOCB(L,K)=0.0  
            ENDDO  
          ENDDO  
          IF(IREAD.EQ.0)THEN  
            IF(ISALTYP.EQ.0)THEN  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO) STPOCB(L,1)  
                IF(ISO.GT.0) GOTO 856  
                DO K=2,KB  
                  STPOCB(L,K)=STPOCB(L,1)  
                ENDDO  
              ENDDO  
            ELSE  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,STPOCB(L,1)  
                IF(ISO.GT.0) GOTO 856  
                DO K=2,KB  
                  STPOCB(L,K)=STPOCB(L,1)  
                ENDDO  
              ENDDO  
            ENDIF  
          ENDIF  
          IF(IREAD.EQ.1)THEN  
            IF(ISALTYP.EQ.0)THEN  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO) (STPOCB(L,K),K=1,KB)  
                IF(ISO.GT.0) GOTO 856  
              ENDDO  
            ELSE  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &              (STPOCB(L,K),K=1,KB)  
                IF(ISO.GT.0) GOTO 856  
              ENDDO  
            ENDIF  
          ENDIF  
          IF(IREAD.EQ.2)THEN  
            IF(ISALTYP.EQ.0)THEN  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO) (STPOCB(L,K),K=1,KBINPUT)  
                IF(ISO.GT.0) GOTO 856  
                DO K=KBINPUT,KB  
                  STPOCB(L,K)=STPOCB(L,KBINPUT)  
                ENDDO  
              ENDDO  
            ELSE  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &              (STPOCB(L,K),K=1,KBINPUT)  
                IF(ISO.GT.0) GOTO 856  
                DO K=KBINPUT,KB  
                  STPOCB(L,K)=STPOCB(L,KBINPUT)  
                ENDDO  
              ENDDO  
            ENDIF  
          ENDIF  
          CLOSE(1)  
        ENDIF  
      ENDIF  
C  
C **  OPEN FILE FPOCB.INP TO READ SPATIALY VARYING, TIME CONSTANT  
C **  PARTICULATE ORGANIC CARBON FRACTION FOR EACH SEDIMENT CLASS  
C **  IN BED  
C  
      IVAL=0  
      DO NT=1,NTOX  
        IF(ISTOC(NT).EQ.2.OR.ISTOC(NT).EQ.3)IVAL=1  
      ENDDO  
      IF(IVAL.EQ.1)THEN  
        IF(ISTPOCB.EQ.3)THEN  
          PRINT *,'READING FPOCB.INP'
          OPEN(1,FILE='FPOCB.INP',STATUS='UNKNOWN')  
          DO NS=1,NSED+NSND  
            DO IS=1,8  
              READ(1,*)  
            ENDDO  
            READ(1,*)ISALTYP,IREAD,KBINPUT  
C
            IF(IREAD.EQ.0)THEN  
              IF(ISALTYP.EQ.0)THEN  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO) STFPOCB(L,1,NS)  
                  IF(ISO.GT.0) GOTO 856  
                  DO K=2,KB  
                    STFPOCB(L,K,NS)=STFPOCB(L,1,NS)  
                  ENDDO  
                ENDDO  
              ELSE  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,STFPOCB(L,1,NS)  
                  IF(ISO.GT.0) GOTO 856  
                  DO K=2,KB  
                    STFPOCB(L,K,NS)=STFPOCB(L,1,NS)  
                  ENDDO  
                ENDDO  
              ENDIF  
            ENDIF  
            IF(IREAD.EQ.1)THEN  
              IF(ISALTYP.EQ.0)THEN  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO) (STFPOCB(L,K,NS),K=1,KB)  
                  IF(ISO.GT.0) GOTO 856  
                ENDDO  
              ELSE  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &                (STFPOCB(L,K,NS),K=1,KB)  
                  IF(ISO.GT.0) GOTO 856  
                ENDDO  
              ENDIF  
            ENDIF  
            IF(IREAD.EQ.2)THEN  
              IF(ISALTYP.EQ.0)THEN  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO) (STFPOCB(L,K,NS),K=1,KBINPUT)  
                  IF(ISO.GT.0) GOTO 856  
                  DO K=KBINPUT,KB  
                    STFPOCB(L,K,NS)=STFPOCB(L,KBINPUT,NS)  
                  ENDDO  
                ENDDO  
              ELSE  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &                (STFPOCB(L,K,NS),K=1,KBINPUT,NS)  
                  IF(ISO.GT.0) GOTO 856  
                  DO K=KBINPUT,KB  
                    STFPOCB(L,K,NS)=STFPOCB(L,KBINPUT,NS)  
                  ENDDO  
                ENDDO  
              ENDIF  
            ENDIF  
          ENDDO  
          CLOSE(1)  
        ENDIF  
      ENDIF  
C**********************************************************************C
C###########################################################################
C HQI Change to include sptially varying, but time constant bulk foc 
C FPOCB  - Bulk foc from data
C PFPOCB - Pseudo foc from data, to be used for all partitioning calculations
C RM, 02/29/04
C**********************************************************************C
C
C **  OPEN FILE FOCB.INP TO READ SPATIALY VARYING, TIME CONSTANT 
C **  PARTICULATE ORGANIC CARBON IN BED AND PSEUDO-POC IN BED
C
      IF(ISTPOCB.EQ.4)THEN
C
        PRINT *,'READING FOCB.INP'
        OPEN(1,FILE='FOCB.INP',STATUS='UNKNOWN')
        DO IS=1,8
          READ(1,*)
        ENDDO
        READ(1,*)ISALTYP,IREAD,KBINPUT
        DO K=1,KB
          DO L=2,LC-1
            FPOCB(L,K)=0.0
          ENDDO
        ENDDO
        DO L=2,LC-1
           READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(FPOCB(L,K),K=1,KBINPUT)
           DO K=1,KBINPUT
              FPOCB(L,K) = FPOCB(L,K)/1000000.
           ENDDO
           DO K=KBINPUT+1,KB
              FPOCB(L,K)=FPOCB(L,KBINPUT)
           ENDDO
           IF(ISO.GT.0) GOTO 856
        ENDDO
        CLOSE(1)
C
        PRINT *,'READING PSEUDO_FOCB.INP'
        OPEN(1,FILE='PSEUDO_FOCB.INP',STATUS='UNKNOWN')
        DO IS=1,8
          READ(1,*)
        ENDDO
        READ(1,*)ISALTYP,IREAD,KBINPUT
        DO L=2,LC-1
           READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(PFPOCB(L,K),K=1,KBINPUT)
           DO K=1,KBINPUT
              PFPOCB(L,K) = PFPOCB(L,K)/1000000.
           ENDDO
           DO K=KBINPUT+1,KB
              PFPOCB(L,K)=PFPOCB(L,KBINPUT)
           ENDDO
           IF(ISO.GT.0) GOTO 856
        ENDDO
        CLOSE(1)
      ENDIF
C  
C **  READ IN INITIAL SALINITY, TEMPERATURE, DYE, SED, SND, TOX  
C **  FOR COLD STARTS FORM FILE XXXX.INP  
C **  SALINITY  
C  
      DO K=1,KC  
        DO L=2,LA  
          SALINIT(L,K)=0.  
        ENDDO  
      ENDDO  

      IF(ISTRAN(1).GE.1.AND.(ISRESTI.EQ.0                  .OR.  
     &                      (ISRESTI.GE.1.AND.ISCI(1).EQ.0).OR.
     &                      (ISTOPT(1).GT.1)))THEN  ! *** PMC SINGLE LINE - FORCE IC
        IF(ISTOPT(1).GE.1)THEN  
          IF (PARTID == MASTER_TASK) PRINT *,'READING SALT.INP'
          OPEN(1,FILE='SALT.INP',STATUS='UNKNOWN')  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,4  
            READ(1,*)  
          ENDDO  
          READ(1,*)ISALTYP  
          IF(ISALTYP.EQ.0)THEN  
            DO L=2,LC-1  
              READ(1,*,IOSTAT=ISO) (SALINIT(L,K),K=1,KC)  
              IF(ISO.GT.0) GOTO 840  
            ENDDO  
          ELSE  
            DO L=2,LC_GLOBAL-1  
              READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(CONINIT(K),K=1,KC) 
              I = XLOC(IDUM)
              J = YLOC(JDUM)
              LMAP = LIJ(I,J)
              IF(LMAP >= 2 .AND. LMAP<=LA)THEN
                SALINIT(LMAP,1:KC) = CONINIT(1:KC) 
              ENDIF
              IF(ISO.GT.0) GOTO 840  
            ENDDO  
          ENDIF  
          CLOSE(1)  
        ENDIF  
      ENDIF  
C  
C **  TEMPERATURE  
C  
      DO K=1,KC  
        DO L=2,LA  
          TEMINIT(L,K)=TEMO  
        ENDDO  
      ENDDO  

      IF(ISTRAN(2).GE.1.AND.(ISRESTI.EQ.0                  .OR.  
     &                      (ISRESTI.GE.1.AND.ISCI(2).EQ.0).OR.
     &                      (ISTOPT(2).GT.9)))THEN  ! *** PMC SINGLE LINE - FORCE IC
        IF(ISTOPT(2).GE.1.OR.INITTEMP.GT.0)THEN  
          IF (PARTID == 0) PRINT *,'READING TEMP.INP'
          OPEN(1,FILE='TEMP.INP',STATUS='UNKNOWN')  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,4  
            READ(1,*)  
          ENDDO  
          READ(1,*)ISALTYP  
          IF(ISALTYP.EQ.0)THEN  
            DO L=2,LC-1  
              READ(1,*,IOSTAT=ISO) (TEMINIT(L,K),K=1,KC)  
              IF(ISO.GT.0) GOTO 842  
            ENDDO  
          ELSE  
            DO L=2,LC_GLOBAL-1  
              READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(CONINIT(K),K=1,KC) 
              I = XLOC(IDUM)
              J = YLOC(JDUM)
              IF (I >0 .AND. I <= IC) THEN
                IF (J > 0 .AND. J <= JC )THEN  
                  LMAP = LIJ(I,J)
                  IF(LMAP >= 2 .AND. LMAP<=LA)THEN
                    TEMINIT(LMAP,1:KC) = CONINIT(1:KC) 
                  END IF
                  IF(ISO.GT.0) GOTO 842  
                END IF
              END IF 
            ENDDO  
          ENDIF  
          CLOSE(1)  
        ENDIF  
      ENDIF  
C  
C **  DYE  
C  
      DO K=1,KC  
        DO L=2,LA  
          DYEINIT(L,K)=0.  
        ENDDO  
      ENDDO  
      IF(ISTRAN(3).GE.1.AND.(ISRESTI.EQ.0                  .OR.  
     &                      (ISRESTI.GE.1.AND.ISCI(3).EQ.0).OR.
     &                      (ISTOPT(3).GT.1)))THEN  ! *** PMC SINGLE LINE - FORCE IC
        IF(ISTOPT(3).GE.1)THEN  
          PRINT *,'READING DYE.INP'
          OPEN(1,FILE='DYE.INP',STATUS='UNKNOWN')  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,4  
            READ(1,*)  
          ENDDO  
          READ(1,*)ISALTYP  
          IF(ISALTYP.EQ.0)THEN  
            DO L=2,LC_GLOBAL-1  
              READ(1,*,IOSTAT=ISO) (DYEINIT(L,K),K=1,KC)  
              IF(ISO.GT.0) GOTO 844  
            ENDDO  
          ELSE  
            DO L=2,LC_GLOBAL-1  
              READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(CONINIT(K),K=1,KC) 
              I = XLOC(IDUM)
              J = YLOC(JDUM)
              IF (I > 0 .AND. I <= IC .AND. J > 0 .AND. J <=JC) THEN
                LMAP = LIJ(I,J)
                IF(LMAP >= 2 .AND. LMAP<=LA)THEN
                  DYEINIT(LMAP,1:KC) = CONINIT(1:KC) 
                END IF
              END IF
              IF(ISO.GT.0) GOTO 844  
            ENDDO  
          ENDIF  
          CLOSE(1)  
        ENDIF  
      ENDIF  
C  
C **  SFL  
C  
      DO K=1,KC  
        DO L=2,LA  
          SFLINIT(L,K)=0.  
        ENDDO  
      ENDDO  
      IF(ISRESTI.EQ.0.AND.ISTRAN(4).GE.1)THEN  
        IF(ISTOPT(4).GE.1)THEN  
          PRINT *,'READING SFL.INP'
          OPEN(1,FILE='SFL.INP',STATUS='UNKNOWN')  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,4  
            READ(1,*)  
          ENDDO  
          READ(1,*)ISALTYP  
          IF(ISALTYP.EQ.0)THEN  
            DO L=2,LC-1  
              READ(1,*,IOSTAT=ISO) (SFLINIT(L,K),K=1,KC)  
              IF(ISO.GT.0) GOTO 846  
            ENDDO  
          ELSE  
            DO L=2,LC-1  
              READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(SFLINIT(L,K),K=1,KC)  
              IF(ISO.GT.0) GOTO 846  
            ENDDO  
          ENDIF  
          CLOSE(1)  
        ENDIF  
      ENDIF  
C  
C **  TOXICS  
C  
      IF(ISTRAN(5).EQ.0) THEN  
        DO NT=1,NTOX  
          DO K=1,KC  
            DO L=2,LA  
              TOXINIT(L,K,NT)=0.0  
            ENDDO  
          ENDDO  
        ENDDO  
        DO NT=1,NTOX  
          DO K=1,KB  
            DO L=2,LA  
              TOXBINIT(L,K,NT)=0.0  
            ENDDO  
          ENDDO  
        ENDDO  
      ENDIF  
      IF(ISTRAN(5).GE.1) THEN  
        DO NT=1,NTOX  
          DO K=1,KC  
            DO L=2,LA  
              TOXINIT(L,K,NT)=TOXINTW(NT)  
            ENDDO  
          ENDDO  
        ENDDO  
        DO NT=1,NTOX  
          DO K=1,KB  
            DO L=2,LA  
              TOXBINIT(L,K,NT)=TOXINTB(NT)  
            ENDDO  
          ENDDO  
        ENDDO  
        IISTMP=1  
        IF(ISRESTI.EQ.0) IISTMP=0  
        IF(ISRESTI.GE.1.AND.ISCI(5).EQ.0) IISTMP=0  
        IF(IISTMP.EQ.0.AND.ISTRAN(5).GE.1)THEN  
          IF(ISLTMT.EQ.0)THEN  
            PRINT *,'READING TOXW.INP'
            OPEN(1,FILE='TOXW.INP',STATUS='UNKNOWN')  
            IF(ITXINT(1).EQ.1.OR.ITXINT(1).EQ.3)THEN  
              DO NT=1,NTOX  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
                DO IS=1,8  
                  READ(1,*)  
                ENDDO  
                READ(1,*)ISALTYP,ITOXWU(NT)  
                IF(ISALTYP.EQ.0)THEN  
                  DO L=2,LC-1  
                    READ(1,*,IOSTAT=ISO) (TOXINIT(L,K,NT),K=1,KC)  
                    IF(ISO.GT.0) GOTO 848  
                  ENDDO  
                ELSE  
                  DO L=2,LC-1  
                    READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(TOXINIT(L,K,
     &                  NT),K=1,KC)  
                    IF(ISO.GT.0) GOTO 848  
                  ENDDO  
                ENDIF  
              ENDDO  
            ENDIF  
            CLOSE(1)  
          ENDIF  
        ENDIF  
        IISTMP=1  
        IF(ISRESTI.EQ.0) IISTMP=0  
        IF(ISRESTI.GE.1.AND.ISCI(5).EQ.0) IISTMP=0  
        IF(IISTMP.EQ.0.AND.ISTRAN(5).GE.1)THEN  
          IF(ISLTMT.EQ.0.)THEN  
            PRINT *,'READING TOXB.INP'
            OPEN(1,FILE='TOXB.INP',STATUS='UNKNOWN')  
            IF(ITXINT(1).EQ.2.OR.ITXINT(1).EQ.3)THEN  
              DO NT=1,NTOX  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
                DO IS=1,8  
                  READ(1,*)  
                ENDDO  
                READ(1,*)ISALTYP,ITOXBU(NT),KBINPUT  
                IF(ISALTYP.EQ.0)THEN  
                  DO L=2,LC-1  
                    READ(1,*,IOSTAT=ISO) (TOXBINIT(L,K,NT),K=1,KBINPUT)  
                    IF(ISO.GT.0) GOTO 852  
                  ENDDO  
                ELSE  
                  DO L=2,LC-1  
                    READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &                  (TOXBINIT(L,K,NT),K=1,KBINPUT)  
                    IF(ISO.GT.0) GOTO 852  
                  ENDDO  
                ENDIF  
              ENDDO  
            ENDIF  
            CLOSE(1)  
          ENDIF  
        ENDIF  
      ENDIF  
C  
C **  COHESIVE SEDIMENT  
C  
      IF(ISTRAN(6).EQ.0)THEN  
        DO NS=1,NSED  
          DO K=1,KC  
            DO L=2,LA  
              SEDINIT(L,K,NS)=0.0  
            ENDDO  
          ENDDO  
        ENDDO  
        DO NS=1,NSED  
          DO K=1,KB  
            DO L=2,LA  
              SEDBINIT(L,K,NS)=0.0  
            ENDDO  
          ENDDO  
        ENDDO  
      END IF  
      IF(ISTRAN(6).GE.1)THEN  
        DO NS=1,NSED  
          DO K=1,KC  
            DO L=2,LA  
              SEDINIT(L,K,NS)=SEDO(NS)  
            ENDDO  
          ENDDO  
        ENDDO  
        DO NS=1,NSED  
          DO K=1,KB  
            DO L=2,LA  
              SEDBINIT(L,K,NS)=SEDBO(NS)  
            ENDDO  
          ENDDO  
        ENDDO  
        ITXINTT=0  
        IF(ISEDINT.EQ.1) ITXINTT=1  
        IF(ISEDINT.EQ.3) ITXINTT=1  
        IISTMP=1  
        IF(ISRESTI.EQ.0) IISTMP=0  
        IF(ISRESTI.GE.1.AND.ISCI(6).EQ.0) IISTMP=0  
        IF(IISTMP.EQ.0.AND.ISTRAN(6).GE.1)THEN  
          IF(ITXINTT.GE.1)THEN  
            PRINT *,'READING SEDW.INP'
            OPEN(1,FILE='SEDW.INP',STATUS='UNKNOWN')  
            DO NS=1,NSED  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
              DO IS=1,8  
                READ(1,*)  
              ENDDO  
              READ(1,*)ISALTYP  
              IF(ISALTYP.EQ.0)THEN  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO) (SEDINIT(L,K,NS),K=1,KC)  
                  IF(ISO.GT.0) GOTO 854  
                ENDDO  
              ELSE  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(SEDINIT(L,K,NS)
     &                ,K=1,KC)  
                  IF(ISO.GT.0) GOTO 854  
                ENDDO  
              ENDIF  
            ENDDO  
            CLOSE(1)  
          ENDIF  
        ENDIF  
        ITXINTT=0  
        IF(ISEDINT.EQ.2) ITXINTT=1  
        IF(ISEDINT.EQ.3) ITXINTT=1  
        IISTMP=1  
        IF(ISRESTI.EQ.0) IISTMP=0  
        IF(ISRESTI.GE.1.AND.ISCI(6).EQ.0) IISTMP=0  
        IF(IISTMP.EQ.0.AND.ISTRAN(6).GE.1.AND.IWRSP(1)/=98)THEN  !avoids loop if SEDZLJ is active
          IF(ITXINTT.GE.1)THEN  
            PRINT *,'READING SEDB.INP'
            OPEN(1,FILE='SEDB.INP',STATUS='UNKNOWN')  
            DO NS=1,NSED  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
              DO IS=1,8  
                READ(1,*)  
              ENDDO  
              READ(1,*)ISALTYP,ISEDBU(NS),KBINPUT  
              DO K=1,KB  
                DO L=2,LC-1  
                  SEDBINIT(L,K,NS)=0.0  
                ENDDO  
              ENDDO  
              IF(ISALTYP.EQ.0)THEN  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO) (SEDBINIT(L,K,NS),K=1,KBINPUT)  
                  IF(ISO.GT.0) GOTO 856  
                ENDDO  
              ELSE  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &                (SEDBINIT(L,K,NS),K=1,KBINPUT)  
                  IF(ISO.GT.0) GOTO 856  
                ENDDO  
              ENDIF  
            ENDDO  
            CLOSE(1)  
          ENDIF  
        ENDIF  
      ENDIF  
C  
C **  NON-COHESIVE SEDIMENT  
C  
      IF(ISTRAN(7).EQ.0)THEN  
        DO NX=1,NSND  
          DO K=1,KC  
            DO L=2,LA  
              SNDINIT(L,K,NX)=0.0  
            ENDDO  
          ENDDO  
        ENDDO  
        DO NX=1,NSND  
          DO K=1,KB  
            DO L=2,LA  
              SNDBINIT(L,K,NX)=0.0  
            ENDDO  
          ENDDO  
        ENDDO  
      END IF  
      IF(ISTRAN(7).GE.1)THEN  
        DO NX=1,NSND  
          NS=NX+NSED  
          DO K=1,KC  
            DO L=2,LA  
              SNDINIT(L,K,NX)=SEDO(NS)  
            ENDDO  
          ENDDO  
        ENDDO  
        DO NX=1,NSND  
          NS=NX+NSED  
          DO K=1,KB  
            DO L=2,LA  
              SNDBINIT(L,K,NX)=SEDBO(NS)  
            ENDDO  
          ENDDO  
        ENDDO  
        ITXINTT=0  
        IF(ISEDINT.EQ.1) ITXINTT=1  
        IF(ISEDINT.EQ.3) ITXINTT=1  
        IISTMP=1  
        IF(ISRESTI.EQ.0) IISTMP=0  
        IF(ISRESTI.GE.1.AND.ISCI(7).EQ.0) IISTMP=0  
        IF(IISTMP.EQ.0.AND.ISTRAN(7).GE.1)THEN  
          IF(ITXINTT.GE.1)THEN  
            PRINT *,'READING SNDW.INP'
            OPEN(1,FILE='SNDW.INP',STATUS='UNKNOWN')  
            DO NX=1,NSND  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
              DO IS=1,8  
                READ(1,*)  
              ENDDO  
              READ(1,*)ISALTYP  
              IF(ISALTYP.EQ.0)THEN  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO) (SNDINIT(L,K,NX),K=1,KC)  
                  IF(ISO.GT.0) GOTO 858  
                ENDDO  
              ELSE  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,(SNDINIT(L,K,NX)
     &                ,K=1,KC)  
                  IF(ISO.GT.0) GOTO 858  
                ENDDO  
              ENDIF  
            ENDDO  
            CLOSE(1)  
          ENDIF  
        ENDIF  
        ITXINTT=0  
        IF(ISEDINT.EQ.2) ITXINTT=1  
        IF(ISEDINT.EQ.3) ITXINTT=1  
        IISTMP=1  
        IF(ISRESTI.EQ.0) IISTMP=0  
        IF(ISRESTI.GE.1.AND.ISCI(7).EQ.0) IISTMP=0  
        IF(IISTMP.EQ.0.AND.ISTRAN(7).GE.1)THEN  
          IF(ITXINTT.GE.1)THEN  
            PRINT *,'READING SNDB.INP'
            OPEN(1,FILE='SNDB.INP',STATUS='UNKNOWN')  
            DO NX=1,NSND  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
              DO IS=1,8  
                READ(1,*)  
              ENDDO  
              DO K=1,KB  
                DO L=2,LC-1  
                  SNDBINIT(L,K,NX)=0.0  
                ENDDO  
              ENDDO  
              READ(1,*)ISALTYP,ISNDBU(NX),KBINPUT  
              IF(ISALTYP.EQ.0)THEN  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO) (SNDBINIT(L,K,NX),K=1,KBINPUT)  
                  IF(ISO.GT.0) GOTO 862  
                ENDDO  
              ELSE  
                DO L=2,LC-1  
                  READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &                (SNDBINIT(L,K,NX),K=1,KBINPUT)  
                  IF(ISO.GT.0) GOTO 862  
                ENDDO  
              ENDIF  
            ENDDO  
            CLOSE(1)  
          ENDIF  
        ENDIF  
      ENDIF  
C  
C  ** SEDIMENT BED MECHANICAL INITIAL CONDITIONS  
C  
      IISTMP=1  
      IF(ISRESTI.EQ.0) IISTMP=0  
      IF(ISRESTI.GE.1.AND.ISCI(6).EQ.0) IISTMP=0  
      IF(ISRESTI.GE.1.AND.ISCI(7).EQ.0) IISTMP=0  
      IF(ISTRAN(6).GT.0.OR.ISTRAN(7).GT.0) THEN  
C  
C  ** BED LAYER THICKNESS  
C  
        IF(IISTMP.EQ.0.AND.IBMECH.GE.1)THEN  
          PRINT *,'READING BEDLAY.INP'
          OPEN(1,FILE='BEDLAY.INP',STATUS='UNKNOWN')  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,8  
            READ(1,*)  
          ENDDO  
          READ(1,*)IBEDLAYU,ISALTYP,KBINPUT  
          IF(IBEDLAYU.GT.0) THEN  
            DO K=1,KB  
              DO L=2,LC-1  
                BEDLINIT(L,K)=0.0  
              ENDDO  
            ENDDO  
            IF(ISALTYP.EQ.0)THEN  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO) (BEDLINIT(L,K),K=1,KBINPUT)  
                IF(ISO.GT.0) GOTO 862  
              ENDDO  
            ELSE  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &              (BEDLINIT(L,K),K=1,KBINPUT)  
                IF(ISO.GT.0) GOTO 862  
              ENDDO  
            ENDIF  
          ENDIF  
          CLOSE(1)  
        ENDIF  
C  
C  ** BED LAYER BULK DENSITY  
C  
        IF(IISTMP.EQ.0.AND.IBMECH.GE.1)THEN  
          PRINT *,'READING BEDBDN.INP'
          OPEN(1,FILE='BEDBDN.INP',STATUS='UNKNOWN')  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,8  
            READ(1,*)  
          ENDDO  
          READ(1,*)IBEDBDNU,ISALTYP,KBINPUT  
          IF(IBEDBDNU.GT.0)THEN  
            DO K=1,KB  
              DO L=2,LC-1  
                BEDBINIT(L,K)=0.0  
              ENDDO  
            ENDDO  
            IF(ISALTYP.EQ.0)THEN  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO) (BEDBINIT(L,K),K=1,KBINPUT)  
                IF(ISO.GT.0) GOTO 862  
              ENDDO  
            ELSE  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &              (BEDBINIT(L,K),K=1,KBINPUT)  
                IF(ISO.GT.0) GOTO 862  
              ENDDO  
            ENDIF  
          ENDIF  
          CLOSE(1)  
        ENDIF  
C  
C  ** BED LAYER DRY DENSITY, POROSITY OR VOID RATIO  
C  
        IF(IISTMP.EQ.0.AND.IBMECH.GE.1)THEN  
          PRINT *,'READING BEDDDN.INP'
          OPEN(1,FILE='BEDDDN.INP',STATUS='UNKNOWN')  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,8  
            READ(1,*)  
          ENDDO  
          READ(1,*)IBEDDDNU,ISALTYP,KBINPUT  
          IF(IBEDDDNU.GT.0)THEN  
            DO K=1,KB  
              DO L=2,LC-1  
                BEDDINIT(L,K)=0.0  
              ENDDO  
            ENDDO  
            IF(ISALTYP.EQ.0)THEN  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO) (BEDDINIT(L,K),K=1,KBINPUT)  
                IF(ISO.GT.0) GOTO 862  
              ENDDO  
            ELSE  
              DO L=2,LC-1  
                READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,  
     &              (BEDDINIT(L,K),K=1,KBINPUT)  
                IF(ISO.GT.0) GOTO 862  
              ENDDO  
            ENDIF  
          ENDIF  
          CLOSE(1)  
        ENDIF  
C  
C  ** CONSOLIDATION MAP  
C  
        IF(IBMECH.EQ.9)THEN  
          PRINT *,'READING CONSOLMAP.INP'
          OPEN(1,FILE='CONSOLMAP.INP',STATUS='UNKNOWN')  
C  
C **    SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,6  
            READ(1,*)  
          ENDDO  
          READ(1,*)ISALTYP  
          IF(ISALTYP.EQ.0)THEN  
            DO L=2,LC-1  
              READ(1,*,IOSTAT=ISO) LCONSOL(L)  
              IF(ISO.GT.0) GOTO 862  
            ENDDO  
          ELSE  
            DO L=2,LC-1  
              READ(1,*,IOSTAT=ISO)LDUM,IDUM,JDUM,LCONSOL(L)  
              IF(ISO.GT.0) GOTO 862  
            ENDDO  
          ENDIF  
          CLOSE(1)  
        ENDIF  
      ENDIF  
   19 FORMAT (/,' INITIAL BUOYANCY ARRAY:',//)  
  907 FORMAT(12F6.2)
      
**********************************************************************C
C
C **  READ IN HORIZONTAL VELOCITY TIME SERIES FOR DATA ASSIMILATION
C **  FILE UVSER.INP
C
C----------------------------------------------------------------------C
C
      IF(NUVSER.GE.1)THEN
        INQUIRE(FILE='UVSER.INP', EXIST=EXIST)
        IF(EXIST)THEN
          OPEN(1,FILE='UVSER.INP',STATUS='UNKNOWN')
        ELSE
          PRINT*,"NUVSER>=1, BUT NO UVSER.INP FILE"
          PAUSE
          STOP
        ENDIF
      
C
C **  SKIP OVER TITLE AND AND HEADER LINES
C 
        DO IS=1,19
          READ(1,'(A120)')
        ENDDO
6669  FORMAT('READING UVSER.INP: NS,MUVSER,TC,TA,RM,AD',1X,I4,1X,I2,1X,4F9.2)
C
        DO NS=1,NUVSER
          READ(1,*,IOSTAT=ISO) ITYPE,MUVSER(NS),TCUVSER(NS),TAUVSER(NS),
     &                   RMULADJ,ADDADJ
	    IF(PARTID==MASTER_TASK)WRITE(6,6669)NS,MUVSER(NS),TCUVSER(NS),
     &                   TAUVSER(NS),RMULADJ,ADDADJ
          IF(ISO.GT.0) GOTO 850
          IF(ITYPE.EQ.-1)THEN
            READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)
            IF(ISO.GT.0) GOTO 860
            DO M=1,MUVSER(NS)
              READ(1,*,IOSTAT=ISO)TUVSER(M,NS),USERTMP,VSERTMP
              IF(ISO.GT.0) GOTO 860
              TUVSER(M,NS)=TUVSER(M,NS)+TAUVSER(NS)
              USERTMP=(RMULADJ*(USERTMP+ADDADJ))
              VSERTMP=(RMULADJ*(VSERTMP+ADDADJ))
	        USER(M,1:KC,NS)=USERTMP*WKQ(1:KC)
		    VSER(M,1:KC,NS)=VSERTMP*WKQ(1:KC)
            ENDDO
          ELSEIF(ITYPE.EQ.0)THEN
            DO M=1,MUVSER(NS)
              READ(1,*,IOSTAT=ISO)TUVSER(M,NS),(USER(M,K,NS),K=1,KC),
     &                                      (VSER(M,K,NS),K=1,KC)
              IF(ISO.GT.0) GOTO 860
              TUVSER(M,NS)=TUVSER(M,NS)+TAUVSER(NS)
              USER(M,1:KC,NS)=RMULADJ*(USER(M,1:KC,NS)+ADDADJ)
              VSER(M,1:KC,NS)=RMULADJ*(VSER(M,1:KC,NS)+ADDADJ)
            ENDDO
          ELSEIF(ITYPE.GE.1)THEN
            DO M=1,MUVSER(NS)
              READ(1,*,IOSTAT=ISO)TUVSER(M,NS),USERTMP,VSERTMP
              IF(ISO.GT.0) GOTO 860
              TUVSER(M,NS)=TUVSER(M,NS)+TAUVSER(NS)
              USERTMP=(RMULADJ*(USERTMP+ADDADJ))
              VSERTMP=(RMULADJ*(VSERTMP+ADDADJ))
              USER(M,1:KC,NS)=-99999.9
              VSER(M,1:KC,NS)=-99999.9
	        USER(M,ITYPE,NS)=USERTMP
              VSER(M,ITYPE,NS)=VSERTMP
            ENDDO
          ENDIF
        ENDDO
C
        CLOSE(1)

C----------------------------------------------------------------------C
CGR 5/10/2011
 !     print*,'ISACDA set to 1';ISACDA=1;pause    
        IF(ISACDA.EQ.2)THEN
          ALLOCATE(IDC(LCM),JDC(LCM),CDAWT(LCM))
          OPEN(1,FILE='CDAWTS.INP')
          DO IS=1,4
            READ(1,'(A120)')
          ENDDO
          DO L=2,LA
            READ(1,*)LDC,IDC(L),JDC(L),CDAWT(L)
  		  CDAWT(L)=TSCDA*CDAWT(L)
          ENDDO
          CLOSE(1)
        ENDIF
CGR 5/10/2011
        IF(ISACDA.EQ.3)THEN
          ALLOCATE(IDC(NLCDA),JDC(NLCDA),CDAWT(NLCDA),CDAWTNS(NLCDA))
          OPEN(1,FILE='CDAWTV.INP')
          DO IS=1,4
            READ(1,'(A120)')
          ENDDO
          DO L=1,NLCDA
            READ(1,*)LDC,IDC(L),JDC(L),CDAWTNS(L),CDAWT(L)
		  CDAWT(L)=TSCDA*CDAWT(L)
          ENDDO
          CLOSE(1)
        ENDIF
      ENDIF
C
C----------------------------------------------------------------------C
C          
C
 7001 FORMAT(2I6,F12.5,50F10.5)
C
C**********************************************************************C
C
C **  READ IN HORIZONTAL VELOCITY SPATIAL FIELD FOR DATA ASSIMILATION
C **  FILE UVDAMAP.INP
C
C----------------------------------------------------------------------C
C
      IF(NLUVDA.GT.0.AND.NUVSER.eq.0)THEN
        OPEN(1,FILE='UVDAMAP.INP',STATUS='UNKNOWN')
        OPEN(99,FILE='UVDAMAP.OUT',STATUS='UNKNOWN')
C
C **  SKIP OVER TITLE AND AND HEADER LINES
C 
        DO IS=1,7
          READ(1,'(A120)')
        ENDDO
C
        UVMAX=0.
	  IMAX=0
	  JMAX=0
        NL=0
        USERT=0.0
        VSERT=0.0
        DO NLT=1,NLUVDA
	    READ(1,*)ITMP,JTMP,USTMP,VSTMP
          ITMP=XLOC(ITMP)
          JTMP=YLOC(JTMP)
          L=LIJ(ITMP,JTMP)
          IF (CELL_INSIDE_DOMAIN(L)) THEN
	      IF(IJCT(ITMP,JTMP).NE.0.AND.IJCT(ITMP,JTMP).NE.9)THEN
	        NL=NL+1
	        ICUVDA(NL)=ITMP
		    JCUVDA(NL)=JTMP
		    USERT(1,NL)=USTMP
		    VSERT(1,NL)=VSTMP
	        UVTMP=SQRT(USTMP*USTMP+VSTMP*VSTMP)
	        WRITE(99,1999)ITMP,JTMP,USTMP,VSTMP,UVTMP
	        IF(UVTMP.GT.UVMAX)THEN
	          UVMAX=UVTMP
	          IMAX=ITMP
	          JMAX=JTMP
	        ENDIF
	        DO K=1,KC
	          USERT(K,NL)=USERT(1,NL)
	          VSERT(K,NL)=VSERT(1,NL)
              ENDDO
              NUVSERA(NL)=NL
              TSUVDA(NL)=TSUVDA(1)
              FSUVDA(NL)=FSUVDA(1)
              IWUVDA(NL)=IWUVDA(1)
              IRVUDA(NL)=IRVUDA(1)
              RRUVDA(NL)=RRUVDA(1)
            ENDIF
          ENDIF
        ENDDO
	  NLUVDA=NL      
	  WRITE(8,*)'I,J,UVMAX=',IMAX,JMAX,UVMAX  
C
        CLOSE(1)
	  CLOSE(99)
      ENDIF
C
 1999 FORMAT(2I6,3F10.3)
C
C**********************************************************************C      
      
C  
C **  READ IN OPEN BOUNDARY SURFACE ELEVATION TIME SERIES FROM THE  
C **  FILE PSER.INP  
C  
      IF(NPSER.GE.1)THEN  
        IF (PARTID ==0) PRINT *,'READING PSER.INP'
        OPEN(1,FILE='PSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,6  
          READ(1,*)  
        ENDDO  
        DO NS=1,NPSER  
          READ(1,*,IOSTAT=ISO) MPSER(NS),TCPSER(NS),TAPSER(NS),  
     &        RMULADJ,ADDADJ  
          IF(ISO.GT.0) GOTO 850  
          DO M=1,MPSER(NS)  
            READ(1,*,IOSTAT=ISO)TPSER(M,NS),PSERTMP  
            IF(ISO.GT.0) GOTO 850  
            TPSER(M,NS)=TPSER(M,NS)+TAPSER(NS)  
            PSER(M,NS)=G*(PSERTMP+ADDADJ)*RMULADJ  
          ENDDO  
        ENDDO  
        CLOSE(1)  
      ENDIF  
 6776 FORMAT(A20)  
C  
C **  READ IN VOLUMETRIC SOURCE OR RIVER INFLOW TIME SERIES FROM THE  
C **  FILE QSER.INP  
C  
      IF(NQSER.GE.1)THEN  
        IF (PARTID ==0 )   PRINT *,'READING QSER.INP'
        OPEN(1,FILE='QSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,14  
          READ(1,*)  
        ENDDO  
        DO NS=1,NQSER  
          READ(1,*,IOSTAT=ISO)ISTYP, MQSER(NS),TCQSER(NS),TAQSER(NS),  
     &        RMULADJ,ADDADJ,ICHGQS  
          IF(ISO.GT.0) GOTO 860  
          IF(ISTYP.EQ.1)THEN  
            READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)  
            IF(ISO.GT.0) GOTO 860  
            DO M=1,MQSER(NS)  
              READ(1,*,IOSTAT=ISO)TQSER(M,NS),QSERTMP  
              IF(ISO.GT.0) GOTO 860  
              TQSER(M,NS)=TQSER(M,NS)+TAQSER(NS)  
              QSERTMP=(RMULADJ*(QSERTMP+ADDADJ))  
              IF(ICHGQS.EQ.1) QSERTMP=MAX(QSERTMP,0.0)  
              IF(ICHGQS.EQ.-1) QSERTMP=MIN(QSERTMP,0.0)  
              DO K=1,KC  
                QSER(M,K,NS)=QSERTMP*WKQ(K)  
              ENDDO  
            ENDDO  
          ELSE  
            DO M=1,MQSER(NS)  
              READ(1,*,IOSTAT=ISO)TQSER(M,NS),(QSER(M,K,NS), K=1,KC)  
              IF(ISO.GT.0) GOTO 860  
              TQSER(M,NS)=TQSER(M,NS)+TAQSER(NS)  
              DO K=1,KC  
                QSER(M,K,NS)=RMULADJ*(QSER(M,K,NS)+ADDADJ)  
                IF(ICHGQS.EQ.1) QSER(M,K,NS)=MAX(QSER(M,K,NS),0.0)  
                IF(ICHGQS.EQ.-1) QSER(M,K,NS)=MIN(QSER(M,K,NS),0.0)  
              ENDDO  
            ENDDO  
          ENDIF  
        ENDDO  
        CLOSE(1)  
      ENDIF  
 2222 FORMAT(2I5,F12.7,F12.4)  

   
C **  READ IN FLOW WITHDRAWL-RETURN FLOW AND CONCENTRATION RISE  
C **  TIME SERIES FROM THE FILE QWRS.INP  
C  
      IF(NQWRSR.GE.1)THEN  
        PRINT *,'READING QWRS.INP'
        OPEN(1,FILE='QWRS.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        NCTMP=4+NSED+NSND+NTOX
        ipmc = NWQV
        IF(ISTRAN(8).GT.0)NCTMP=NCTMP+ipmc   ! *** IF NQWV CHANGES THIS SHOULD BE UPDATED
        DO IS=1,16  
          READ(1,*)  
        ENDDO  
        DO NS=1,NQWRSR  
          READ(1,*,IOSTAT=ISO)ISTYP,MQWRSR(NS),TCQWRSR(NS),TAQWRSR(NS),  
     &        RMULADJ,ADDADJ  
          IF(ISO.GT.0) GOTO 865  
          IF(ISTYP.EQ.0)THEN  
            DO NC=1,NCTMP  
              DO M=1,MQWRSR(NS)  
                CQWRSER(M,NS,NC)=0.  
              ENDDO  
            ENDDO  
            DO M=1,MQWRSR(NS)  
              READ(1,*,IOSTAT=ISO)TQWRSER(M,NS),QWRSER(M,NS)  
              IF(ISO.GT.0) GOTO 865  
              TQWRSER(M,NS)=TQWRSER(M,NS)+TAQWRSR(NS)  
              QWRSER(M,NS)=(RMULADJ*(QWRSER(M,NS)+ADDADJ))  
            ENDDO  
          ELSE  
            DO M=1,MQWRSR(NS)  
              READ(1,*,IOSTAT=ISO)TQWRSER(M,NS),QWRSER(M,NS),  
     &            (CQWRSER(M,NS,NC),NC=1,NCTMP)  
              IF(ISO.GT.0) GOTO 865  
              TQWRSER(M,NS)=TQWRSER(M,NS)+TAQWRSR(NS)  
              QWRSER(M,NS)=(RMULADJ*(QWRSER(M,NS)+ADDADJ))  
            ENDDO  
          ENDIF  
        ENDDO  
        CLOSE(1)  
      ENDIF  
C  
C **  READ IN GROUNDWATER INFLOW/OUTFLOW AND CONCENTRATION TIME  
C **  SERIES FROM THE FILE GWSER.INP  
C  
      IF(ISGWIT.EQ.2)THEN  
        PRINT *,'READING GWSER.INP'
        OPEN(1,FILE='GWSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        NCTMP=4+NSED+NSND+NTOX  
        DO IS=1,14  
          READ(1,*)  
        ENDDO  
        READ(1,*)NGWSER  
        IF(NGWSER.GT.0)THEN  
          DO NS=1,NGWSER  
            READ(1,*,IOSTAT=ISO)MGWSER(NS),TCGWSER(NS),TAGWSER(NS),  
     &          RMULADJ,ADDADJ  
            IF(ISO.GT.0) GOTO 865  
            DO M=1,MGWSER(NS)  
              READ(1,*,IOSTAT=ISO)TGWSER(M,NS),GWSER(M,NS),  
     &            (GWCSER(M,NS,NC),NC=1,NCTMP)  
              IF(ISO.GT.0) GOTO 865  
              TGWSER(M,NS)=TGWSER(M,NS)+TAGWSER(NS)  
              GWSER(M,NS)=(RMULADJ*(GWSER(M,NS)+ADDADJ))  
            ENDDO  
          ENDDO  
        ENDIF  
        CLOSE(1)  
      ENDIF  
C  
C **  READ IN SPATIAL MAPS AND TIME SERIES FOR EXTERNAL SPECIFICATION OF  
C **  PARTICULATE ORGANIC CARBON FOR USE IN TOXIC CONTAMINANT SORPTION  
C **  DISSOLVED ORGANIC CARBON  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C        READ(1,*)  
C        READ(1,*)NOCSER  
C            READ(1,*,IOSTAT=ISO)MOCSER(NS),TCOCSER(NS),TAOCSER(NS),  
C              READ(1,*,IOSTAT=ISO)TOCSER(M,NS),DOCWSER(M,NS),  
C **  READ IN OPEN BOUNDARY OR VOLUMETRIC SOURCE SALINITY TIME SERIES  
C **  FROM THE FILE SSER.INP  
C  
 8888 FORMAT(3I5,2F10.2)  
      IF(NCSER(1).GE.1)THEN  
        PRINT *,'READING SSER.INP'
        OPEN(1,FILE='SSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,15  
          READ(1,*)  
        ENDDO  
        NC=1  
        DO NS=1,NCSER(NC)  
          READ(1,*,IOSTAT=ISO)ISTYP,MCSER(NS,NC),TCCSER(NS,NC),  
     &        TACSER(NS,NC),RMULADJ,ADDADJ  
          IF(ISO.GT.0) GOTO 870  
          IF(ISTYP.EQ.1)THEN  
            READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)  
            IF(ISO.GT.0) GOTO 870  
            DO M=1,MCSER(NS,NC)  
              READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),CSERTMP  
              IF(M.EQ.1)WRITE(8,8888)NC,NS,M,TCSER(M,NS,NC),CSERTMP  
              IF(M.EQ.MCSER(NS,NC))WRITE(8,8888)NC,NS,M,  
     &            TCSER(M,NS,NC),CSERTMP  
              IF(ISO.GT.0) GOTO 870  
              TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
              DO K=1,KC  
                CSER(M,K,NS,NC)=(RMULADJ*(CSERTMP+ADDADJ))*WKQ(K)  
              ENDDO  
            ENDDO  
          ELSE  
            DO M=1,MCSER(NS,NC)  
             READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),(CSER(M,K,NS,NC),K=1,KC)  
              IF(ISO.GT.0) GOTO 870  
              TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
              DO K=1,KC  
                CSER(M,K,NS,NC)=RMULADJ*(CSER(M,K,NS,NC)+ADDADJ)  
              ENDDO  
            ENDDO  
          ENDIF  
        ENDDO  
        CLOSE(1)  
      ENDIF  
C  
C **  READ IN OPEN BOUNDARY OR VOLUMETRIC SOURCE TEMPERATURE TIME  
C **  SERIES FROM THE FILE TSER.INP  
C  


      IF(NCSER(2).GE.1)THEN  
        IF (PARTID ==0) PRINT *,'READING TSER.INP'
        OPEN(1,FILE='TSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,15  
          READ(1,*)  
        ENDDO  
        NC=2  
        DO NS=1,NCSER(NC)  
          READ(1,*,IOSTAT=ISO)ISTYP,MCSER(NS,NC),TCCSER(NS,NC),  
     &        TACSER(NS,NC),RMULADJ,ADDADJ  
         IF(ISO.GT.0) GOTO 880  
           IF(ISTYP.EQ.1)THEN  
             READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)  
             IF(ISO.GT.0) GOTO 880  
             DO M=1,MCSER(NS,NC)  
               READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),CSERTMP  
               IF(M.EQ.1)WRITE(8,8888)NC,NS,M,TCSER(M,NS,NC),CSERTMP  
               IF(M.EQ.MCSER(NS,NC))WRITE(8,8888)NC,NS,M,  
     &            TCSER(M,NS,NC),CSERTMP  
               IF(ISO.GT.0) GOTO 880  
               TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
               DO K=1,KC 
                 CSER(M,K,NS,NC)=(RMULADJ*(CSERTMP+ADDADJ))*WKQ(K)
               ENDDO  
            ENDDO  
          ELSE  
            DO M=1,MCSER(NS,NC)  
             READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),(CSER(M,K,NS,NC),K=1,KC)  
              IF(ISO.GT.0) GOTO 880  
              TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
              DO K=1,KC  
                CSER(M,K,NS,NC)=RMULADJ*(CSER(M,K,NS,NC)+ADDADJ)  
              ENDDO  
            ENDDO  
          ENDIF  
        ENDDO  
        CLOSE(1)  
      ENDIF  
  
C **  READ IN OPEN BOUNDARY OR VOLUMETRIC SOURCE DYE CONCENTRATION  
C **  TIME SERIES FROM THE FILE DSER.INP  
C  
      IF(NCSER(3).GE.1)THEN  
        PRINT *,'READING DSER.INP'
        OPEN(1,FILE='DSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,15  
          READ(1,*)  
        ENDDO  
        NC=3  
        DO NS=1,NCSER(NC)  
          READ(1,*,IOSTAT=ISO)ISTYP,MCSER(NS,NC),TCCSER(NS,NC),  
     &        TACSER(NS,NC),RMULADJ,ADDADJ  
          IF(ISO.GT.0) GOTO 890  
          IF(ISTYP.EQ.1)THEN  
            READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)  
            IF(ISO.GT.0) GOTO 890  
            DO M=1,MCSER(NS,NC)  
              READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),CSERTMP  
              IF(ISO.GT.0) GOTO 890  
              TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
              DO K=1,KC  
                CSER(M,K,NS,NC)=(RMULADJ*(CSERTMP+ADDADJ))*WKQ(K)  
              ENDDO  
            ENDDO  
          ELSE  
            DO M=1,MCSER(NS,NC)  
             READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),(CSER(M,K,NS,NC),K=1,KC)  
              IF(ISO.GT.0) GOTO 890  
              TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
              DO K=1,KC  
                CSER(M,K,NS,NC)=RMULADJ*(CSER(M,K,NS,NC)+ADDADJ)  
              ENDDO  
            ENDDO  
          ENDIF  
        ENDDO  
        CLOSE(1)  
      ENDIF  
C  
C **  READ IN OPEN BOUNDARY OR VOLUMETRIC SOURCE COHESIVE SEDIMENT  
C **  CONCENTRATION TIME SERIES FROM THE FILE SDSER.INP  
C  
      IF(NSED.GT.0)THEN  
        NFSED=MSVSED(1)  
        IF(NCSER(NFSED).GE.1)THEN  
          PRINT *,'READING SDSER.INP'
          OPEN(1,FILE='SDSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,15  
            READ(1,*)  
          ENDDO  
          NC=MSVSED(1)  
          DO NS=1,NCSER(NC)  
            READ(1,*,IOSTAT=ISO)ISTYP,MCSER(NS,NC),TCCSER(NS,NC),  
     &          TACSER(NS,NC),RMULADS(1),ADDADS(1)  
            IF(ISO.GT.0) GOTO 900  
            IF(NSED.GT.1)THEN  
              DO NT=2,NSED  
                READ(1,*,IOSTAT=ISO)RMULADS(NT),ADDADS(NT)  
                IF(ISO.GT.0) GOTO 900  
                NTT=NT-1  
                MCSER(NS,NC+NTT)=MCSER(NS,NC)  
                TCCSER(NS,NC+NTT)=TCCSER(NS,NC)  
                TACSER(NS,NC+NTT)=TACSER(NS,NC)  
              ENDDO  
            ENDIF  
            IF(ISTYP.EQ.1)THEN  
              READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)  
              IF(ISO.GT.0) GOTO 900  
              DO M=1,MCSER(NS,NC)  
                READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),CSERTMP  
                IF(ISO.GT.0) GOTO 900  
                TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
                DO K=1,KC  
                 CSER(M,K,NS,NC)=(RMULADS(1)*(CSERTMP+ADDADS(1)))*WKQ(K)  
                ENDDO  
                DO NT=2,NSED  
                  NTT=NT-1  
                  TCSER(M,NS,NC+NTT)=TCSER(M,NS,NC)  
                  READ(1,*,IOSTAT=ISO)CSERTMP  
                  IF(ISO.GT.0) GOTO 900  
                  DO K=1,KC  
                    CSER(M,K,NS,NC+NTT)  
     &                  =(RMULADS(NT)*(CSERTMP+ADDADS(NT)))*WKQ(K)  
                  ENDDO  
                ENDDO  
              ENDDO  
            ELSE  
              DO M=1,MCSER(NS,NC)  
                READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),(CSER(M,K,NS,NC),
     &              K=1,KC)  
                IF(ISO.GT.0) GOTO 900  
                TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
                DO K=1,KC  
                  CSER(M,K,NS,NC)=RMULADS(1)*(CSER(M,K,NS,NC)+ADDADS(1))  
                ENDDO  
                DO NT=2,NSED  
                  NTT=NT-1  
                  TCSER(M,NS,NC+NTT)=TCSER(M,NS,NC)  
                  READ(1,*,IOSTAT=ISO)(CSER(M,K,NS,NC+NTT), K=1,KC)  
                  IF(ISO.GT.0) GOTO 900  
                  DO K=1,KC  
                    CSER(M,K,NS,NC+NTT)  
     &                  =RMULADS(NT)*(CSER(M,K,NS,NC+NTT)+ADDADS(NT))  
                  ENDDO  
                ENDDO  
              ENDDO  
            ENDIF  
          ENDDO  
          CLOSE(1)  
        ENDIF  
      ENDIF  
C  
C **  CHECK SEDIMENT SERIES  
C **  READ IN OPEN BOUNDARY OR VOLUMETRIC SOURCE NONCOHESIVE SEDIMENT  
C **  CONCENTRATION TIME SERIES FROM THE FILE SNSER.INP  
C  
      IF(NSND.GT.0)THEN  
        NFSND=MSVSND(1)  
        IF(NCSER(NFSND).GE.1)THEN  
          PRINT *,'READING SNSER.INP'
          OPEN(1,FILE='SNSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,15  
            READ(1,*)  
          ENDDO  
          NC=MSVSND(1)  
          DO NS=1,NCSER(NC)  
            READ(1,*,IOSTAT=ISO)ISTYP,MCSER(NS,NC),TCCSER(NS,NC),  
     &          TACSER(NS,NC),RMULADS(1),ADDADS(1)  
            IF(ISO.GT.0) GOTO 902  
            IF(NSND.GT.1)THEN  
              DO NT=2,NSND  
                READ(1,*,IOSTAT=ISO)RMULADS(NT),ADDADS(NT)  
                IF(ISO.GT.0) GOTO 902  
                NTT=NT-1  
                MCSER(NS,NC+NTT)=MCSER(NS,NC)  
                TCCSER(NS,NC+NTT)=TCCSER(NS,NC)  
                TACSER(NS,NC+NTT)=TACSER(NS,NC)  
              ENDDO  
            ENDIF  
            IF(ISTYP.EQ.1)THEN  
              READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)  
              IF(ISO.GT.0) GOTO 902  
              DO M=1,MCSER(NS,NC)  
                READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),CSERTMP  
                IF(ISO.GT.0) GOTO 902  
                TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
                DO K=1,KC  
                 CSER(M,K,NS,NC)=(RMULADS(1)*(CSERTMP+ADDADS(1)))*WKQ(K)  
                ENDDO  
                DO NT=2,NSND  
                  NTT=NT-1  
                  TCSER(M,NS,NC+NTT)=TCSER(M,NS,NC)  
                  READ(1,*,IOSTAT=ISO)CSERTMP  
                  IF(ISO.GT.0) GOTO 902  
                  DO K=1,KC  
                    CSER(M,K,NS,NC+NTT)  
     &                  =(RMULADS(NT)*(CSERTMP+ADDADS(NT)))*WKQ(K)  
                  ENDDO  
                ENDDO  
              ENDDO  
            ELSE  
              DO M=1,MCSER(NS,NC)  
                READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),(CSER(M,K,NS,NC),
     &              K=1,KC)  
                IF(ISO.GT.0) GOTO 902  
                TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
                DO K=1,KC  
                  CSER(M,K,NS,NC)=RMULADS(1)*(CSER(M,K,NS,NC)+ADDADS(1))  
                ENDDO  
                DO NT=2,NSND  
                  NTT=NT-1  
                  TCSER(M,NS,NC+NTT)=TCSER(M,NS,NC)  
                  READ(1,*,IOSTAT=ISO)(CSER(M,K,NS,NC+NTT), K=1,KC)  
                  IF(ISO.GT.0) GOTO 902  
                  DO K=1,KC  
                    CSER(M,K,NS,NC+NTT)  
     &                  =RMULADS(NT)*(CSER(M,K,NS,NC+NTT)+ADDADS(NT))  
                  ENDDO  
                ENDDO  
              ENDDO  
            ENDIF  
          ENDDO  
          CLOSE(1)  
        ENDIF  
      ENDIF  
C  
C **  CHECK SEDIMENT SERIES  
C  
 2001 FORMAT(3I5,2F12.5)  
C  
C **  READ IN OPEN BOUNDARY OR VOLUMETRIC SOURCE TOXIC CONTAMINANT  
C **  CONCENTRATION TIME SERIES FROM THE FILE TXSER.INP  
C  
      IF(NTOX.GT.0)THEN  
        NFTOX=MSVTOX(1)  
        IF(NCSER(NFTOX).GE.1)THEN  
          PRINT *,'READING TXSER.INP'
          OPEN(1,FILE='TXSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
          DO IS=1,15  
            READ(1,*)  
          ENDDO  
          NC=MSVTOX(1)  
          DO NS=1,NCSER(NC)  
            READ(1,*,IOSTAT=ISO)ISTYP,MCSER(NS,NC),TCCSER(NS,NC),  
     &          TACSER(NS,NC),RMULADJ,ADDADJ  
            IF(ISO.GT.0) GOTO 904  
            DO NT=2,NTOX  
              NTT=NT-1  
              MCSER(NS,NC+NTT)=MCSER(NS,NC)  
              TCCSER(NS,NC+NTT)=TCCSER(NS,NC)  
              TACSER(NS,NC+NTT)=TACSER(NS,NC)  
            ENDDO  
            IF(ISTYP.EQ.1)THEN  
              READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)  
              IF(ISO.GT.0) GOTO 904  
              DO M=1,MCSER(NS,NC)  
                READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),CSERTMP  
                IF(ISO.GT.0) GOTO 904  
                TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
                DO K=1,KC  
                  CSER(M,K,NS,NC)=(RMULADJ*(CSERTMP+ADDADJ))*WKQ(K)  
                ENDDO  
                DO NT=2,NTOX  
                  NTT=NT-1  
                  TCSER(M,NS,NC+NTT)=TCSER(M,NS,NC)  
                  READ(1,*,IOSTAT=ISO)CSERTMP  
                  IF(ISO.GT.0) GOTO 904  
                  DO K=1,KC  
                   CSER(M,K,NS,NC+NTT)=(RMULADJ*(CSERTMP+ADDADJ))*WKQ(K)  
                  ENDDO  
                ENDDO  
              ENDDO  
            ELSE  
              DO M=1,MCSER(NS,NC)  
                READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),(CSER(M,K,NS,NC),
     &              K=1,KC)  
                IF(ISO.GT.0) GOTO 904  
                TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
                DO K=1,KC  
                  CSER(M,K,NS,NC)=RMULADJ*(CSER(M,K,NS,NC)+ADDADJ)  
                ENDDO  
                DO NT=2,NTOX  
                  NTT=NT-1  
                  TCSER(M,NS,NC+NTT)=TCSER(M,NS,NC)  
                  READ(1,*,IOSTAT=ISO)(CSER(M,K,NS,NC+NTT), K=1,KC)  
                  IF(ISO.GT.0) GOTO 904  
                  DO K=1,KC  
                    CSER(M,K,NS,NC+NTT)=RMULADJ*(CSER(M,K,NS,NC+NTT)
     &                  +ADDADJ)  
                  ENDDO  
                ENDDO  
              ENDDO  
            ENDIF  
          ENDDO  
          CLOSE(1)  
        ENDIF  
      ENDIF  
C  
C **  READ IN OPEN BOUNDARY OR VOLUMETRIC SOURCE SHELL FISH LARVAE  
C **  TIME SERIES FROM THE FILE SFSER.INP  
C  
      IF(NCSER(4).GE.1)THEN  
        PRINT *,'READING SFSER.INP'
        OPEN(1,FILE='SFSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,15  
          READ(1,*)  
        ENDDO  
        NC=7  
        DO NS=1,NCSER(NC)  
          READ(1,*,IOSTAT=ISO)ISTYP,MCSER(NS,NC),TCCSER(NS,NC),  
     &        TACSER(NS,NC),RMULADJ,ADDADJ  
          IF(ISO.GT.0) GOTO 910  
          IF(ISTYP.EQ.1)THEN  
            READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)  
            IF(ISO.GT.0) GOTO 910  
            DO M=1,MCSER(NS,NC)  
              READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),CSERTMP  
              IF(ISO.GT.0) GOTO 910  
              TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
              DO K=1,KC  
                CSER(M,K,NS,NC)=(RMULADJ*(CSERTMP+ADDADJ))*WKQ(K)  
              ENDDO  
            ENDDO  
          ELSE  
            DO M=1,MCSER(NS,NC)  
             READ(1,*,IOSTAT=ISO)TCSER(M,NS,NC),(CSER(M,K,NS,NC),K=1,KC)  
              IF(ISO.GT.0) GOTO 910  
              TCSER(M,NS,NC)=TCSER(M,NS,NC)+TACSER(NS,NC)  
              DO K=1,KC  
                CSER(M,K,NS,NC)=RMULADJ*(CSER(M,K,NS,NC)+ADDADJ)  
              ENDDO  
            ENDDO  
          ENDIF  
        ENDDO  
        CLOSE(1)  
      ENDIF  
C  
C **  READ IN FREE SURFACE ELEVATION OR PRESSURE CONTROLLED FLOW  
C **  SPECIFICATION FROM THE FILE QCTL.INP  
C **  THE FLOW IS GIVEN BY:  
C         FREE SURFACE  
C         FREE SURFACE  
C        FLOW=0  
C      ELSE  
C        ENTER QCTL(M,K,NS) VS HDIFCTL(M,NS) TABLE WITH DELH TO GIVE  
C  
      IF(NQCTL.GE.1)THEN  
        PRINT *,'READING QCTL.INP'
        OPEN(1,FILE='QCTL.INP',STATUS='UNKNOWN') 
        IF(DEBUG)THEN 
          OPEN(99,FILE='QCTLCK.INP',STATUS='UNKNOWN')  
          CLOSE(99,STATUS='DELETE')  
          OPEN(99,FILE='QCTLCK.INP',STATUS='UNKNOWN')  
        ENDIF
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,14  
          READ(1,*)  
        ENDDO  
        DO NS=1,NQCTLT  
          READ(1,*, IOSTAT=ISO)ISTYP,MQCTL(NS),HCTLUA(NS),HCTLUM(NS),  
     &        HCTLDA(NS),HCTLDM(NS),RMULADJ,ADDADJ,AQCTL(NS)  
          IF(DEBUG)THEN
            WRITE(99,991)NS  
            WRITE(99,992)ISTYP,MQCTL(NS),HCTLUA(NS),HCTLUM(NS),  
     &        HCTLDA(NS),HCTLDM(NS),RMULADJ,ADDADJ,AQCTL(NS)  
          ENDIF
          IF(ISO.GT.0) GOTO 920  
          IF(ISTYP.EQ.0)THEN  
            DO M=1,MQCTL(NS)  
              READ(1,*,IOSTAT=ISO) HDIFCTL(M,NS),(QCTL(M,1,K,NS),K=1,KC)
              IF(ISO.GT.0) GOTO 920  
              DO K=1,KC  
                QCTL(M,1,K,NS)=RMULADJ*(QCTL(M,1,K,NS)+ADDADJ)  
              ENDDO  
            ENDDO  
          ENDIF  
          IF(ISTYP.EQ.1)THEN  
            READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)  
            IF(ISO.GT.0) GOTO 920  
            DO M=1,MQCTL(NS)  
              READ(1,*,IOSTAT=ISO) HDIFCTL(M,NS),QCTLTMP  
              IF(ISO.GT.0) GOTO 920  
              DO K=1,KC  
                QCTL(M,1,K,NS)=RMULADJ*(QCTLTMP+ADDADJ)*WKQ(K)  
              ENDDO  
            ENDDO  
          ENDIF  
          IF(ISTYP.EQ.2)THEN  
            DO MD=1,MQCTL(NS)  
              DO MU=1,MQCTL(NS)  
                READ(1,*,IOSTAT=ISO) HDIFCTL(MU,NS),HDIFCTD(MD,NS),  
     &              (QCTL(MU,MD,K,NS),K=1,KC)  
                IF(ISO.GT.0) GOTO 920  
                DO K=1,KC  
                  QCTL(MU,MD,K,NS)=RMULADJ*(QCTL(MU,MD,K,NS)+ADDADJ)  
                ENDDO  
              ENDDO  
            ENDDO  
          ENDIF  
          IF(ISTYP.EQ.3)THEN  
            READ(1,*,IOSTAT=ISO) (WKQ(K),K=1,KC)  
            IF(ISO.GT.0) GOTO 920  
            DO MD=1,MQCTL(NS)  
              DO MU=1,MQCTL(NS)  
               READ(1,*,IOSTAT=ISO)HDIFCTL(MU,NS),HDIFCTD(MD,NS),QCTLTMP  
                IF(ISO.GT.0) GOTO 920  
                DO K=1,KC  
                  QCTL(MU,MD,K,NS)=RMULADJ*(QCTLTMP+ADDADJ)*WKQ(K)  
                ENDDO  
              ENDDO  
            ENDDO  
          ENDIF
          IF(DEBUG)THEN  
            IF(ISTYP.LE.1)THEN  
              DO M=1,MQCTL(NS)  
                WRITE(99,993)M,HDIFCTL(M,NS),(QCTL(M,1,K,NS),K=1,KC)  
              ENDDO  
            ENDIF  
            IF(ISTYP.GE.2)THEN  
              DO MD=1,MQCTL(NS)  
                DO MU=1,MQCTL(NS)  
                  WRITE(99,994)MU,MD,HDIFCTL(MU,NS),HDIFCTD(MD,NS),  
     &              (QCTL(MU,MD,K,NS),K=1,KC)  
                ENDDO  
              ENDDO  
            ENDIF  
          ENDIF
        ENDDO  
        CLOSE(1)  
        IF(DEBUG)CLOSE(99)  
      ENDIF  
  991 FORMAT(/,'CONTROL TABLE NS =',I5,/)  
  992 FORMAT(2I5,7F10.4)  
  993 FORMAT(I5,11F10.4)  
  994 FORMAT(2I5,11F10.4)  
 1001 FORMAT(/,'READ ERROR FROM FILE EFDC.INP ON CARD ',A3/)  
 1002 FORMAT(/,'INPUT ECHO NCARD = ',A/)  
C  
      DO L=2,LA  
        PATMT(L)=1000.  
        TATMT(L)=0.  
        RAINT(L)=0.  
        EVAPT(L)=0.  
        SOLSWRT(L)=1.  ! *** Address SUNDAY.INP Option
        CLOUDT(L)=0.  
        SVPA(L)=0.  
        RHA(L)=0.  
        VPA(L)=0.  
        CLEVAP(L)=0.  
        CCNHTT(L)=0.  
      ENDDO  
      IF(NASER.GT.0)THEN  
        PRINT *,'READING ASER.INP'
        OPEN(1,FILE='ASER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        ! *** DSLLC BEGIN BLOCK  
        DS_LAT=0.0
        DS_LONG=0.0
        COMPUTESOLRAD=.FALSE.
        DO IS=1,29 
          READ(1,*)  
        ENDDO
        READ(1,22)TEXT
        COMPUTESOLRAD=PARSE_LOGICAL(TEXT)
        READ(1,22)TEXT
        DS_LAT=PARSE_REAL(TEXT)
        READ(1,22)TEXT
        DS_LONG=PARSE_REAL(TEXT)
        READ(1,22)TEXT
        USESHADE=PARSE_LOGICAL(TEXT)

        PRINT *,' NUMBER OF ATMOSPHERE SERIES=',NASER
        PRINT *,' COMPUTESOLRAD=',COMPUTESOLRAD
        PRINT *,' DS_LONG=',DS_LONG
        PRINT *,' DS_LAT=',DS_LAT

        DO IS=1,3
          READ(1,*)  
        ENDDO
        ! *** DSLLC END BLOCK  

        DO N=1,NASER  
          READ(1,*,IOSTAT=ISO) MASER(N),TCASER(N),TAASER(N),IRELH(N),  
     &        RAINCVT,EVAPCVT,SOLRCVT,CLDCVT  
          IF(ISO.GT.0) GOTO 940 
          ! *** These parameters are read in for every series but only the last is
          ! *** actually used 
          READ(1,*,IOSTAT=ISO) IASWRAD,REVC,RCHC,SWRATNF,  
     &        SWRATNS,FSWRATF,DABEDT,TBEDIT,HTBED1,HTBED2  
          IF(ISO.GT.0) GOTO 940  
          DO M=1,MASER(N)  
            READ(1,*,IOSTAT=ISO)TASER(M,N),PATM(M,N),TDRY(M,N),  
     &          TWET(M,N),RAIN(M,N),EVAP(M,N),SOLSWR(M,N),CLOUD(M,N)  
                 IF(ISO.GT.0) GOTO 940
          ENDDO  
          DO M=1,MASER(N)  
            TASER(M,N)=TASER(M,N)+TAASER(N)  
            RAIN(M,N)=RAINCVT*RAIN(M,N)  
            EVAP(M,N)=EVAPCVT*EVAP(M,N)  
            SOLSWR(M,N)=SOLRCVT*SOLSWR(M,N)  
            CLOUD(M,N)=CLDCVT*CLOUD(M,N)
          ENDDO  
        ENDDO  
        CLOSE(1)
        IF(ISTOPT(2).EQ.4.OR.ISTOPT(2).EQ.14)THEN
          ! *** HTBED2 = Bottom Heat Ex Coeff (W / m2 / deg C)
          ! *** RHO    = 1000.0  Density (kg / m^3)  
          ! *** CP     = 4179.0  Specific Heat (J / kg / degC)
          ! *** 0.2393E-6 = 1/RHO/CP
          HTBED2 = HTBED2*0.2393E-6    ! *** m/s

          BETAF=HTBED1
          HTBED1=0.0
        ENDIF  
      ENDIF  
      TEMB(1)=ABS(TBEDIT)
      TEMB(LC)=TEMB(1)
      TEMB1(1)=TEMB(1) 
      TEMB1(LC)=TEMB(1)  
      IF(ISRESTI.EQ.0)THEN  
        DO L=2,LA  
          TEMB(L)=TEMB(1)  
          TEMB1(L)=TEMB(1)  
        ENDDO  
      ENDIF  
C  
C **  READ IN ABOVE WATER SURFACE WIND TIME SERIES FROM THE  
C **  FILE WSER.INP  
C  
      DO L=2,LA  
        WINDST(L)=0.  
        TSX(L)=0.  
        TSY(L)=0.  
      ENDDO  
      IF(NWSER.GT.0)THEN  
        PRINT *,'READING WSER.INP'
        OPEN(1,FILE='WSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,12
          READ(1,*)  
        ENDDO

        READ(1,22)TEXT
        WINDH=PARSE_REAL(TEXT)
        IF(WINDH.LE.0.1)WINDH=10.0

        ! *** PERIOD TO TURN OFF WIND SHEAR
        READ(1,22)TEXT
        WINTER_START=PARSE_REAL(TEXT)
        READ(1,22)TEXT
        WINTER_END=PARSE_REAL(TEXT)
        IF(WINTER_END.LE.WINTER_START)THEN
          WINTER_START = 0.
          WINTER_END=0.
        ENDIF
        IF (PARTID == MASTER_TASK) THEN
          PRINT *,' NUMBER OF WIND SERIES=',NWSER
          PRINT *,' ANEMOMETER HEIGHT (m)=',WINDH
        END IF
        IF (WINTER_START.LT.WINTER_END) THEN
          PRINT *,' SURFACE WIND STRESSES TURNED OFF FROM ',WINTER_START
     &            ," TO ",WINTER_END
        ENDIF
        
        DO IS=1,2
          READ(1,*)  
        ENDDO
  
        DO N=1,NWSER  
          READ(1,*,IOSTAT=ISO) MWSER(N),TCWSER(N),TAWSER(N),WINDSCT,  
     &        ISWDINT(N)  
          IF(ISO.GT.0) GOTO 940  
          DO M=1,MWSER(N)  
            READ(1,*,IOSTAT=ISO)TWSER(M,N),WINDS(M,N),WINDD(M,N)  
            IF(ISO.GT.0) GOTO 940  
          ENDDO  
          DO M=1,MWSER(N)  
            TWSER(M,N)=TWSER(M,N)+TAWSER(N)  
          ENDDO  
          IF(ISWDINT(N).LE.1)THEN  
            DO M=1,MWSER(N)  
              WINDS(M,N)=WINDSCT*WINDS(M,N)  
            ENDDO  
          ENDIF  
          IF(ISWDINT(N).EQ.1)THEN  
            DO M=1,MWSER(N)  
              IF(WINDD(M,N).LE.180.)THEN  
                WINDD(M,N)=WINDD(M,N)+180.  
                IF(WINDD(M,N).EQ.360.) WINDD(M,N)=0.  
              ELSE  
                WINDD(M,N)=WINDD(M,N)-180.  
                IF(WINDD(M,N).EQ.360.) WINDD(M,N)=0.  
              ENDIF  
            ENDDO  
          ENDIF  
          IF(ISWDINT(N).EQ.2)THEN  
            DO M=1,MWSER(N)  
              WINDS(M,N)=WINDSCT*WINDS(M,N)  
              WINDD(M,N)=WINDSCT*WINDD(M,N)  
            ENDDO  
          ENDIF  
        ENDDO  
        CLOSE(1)  
      ENDIF  

C **  READ IN SHELL FISH LARAVE BEHAVIOR DATA  
C **  FROM THE FILE SFBSER.INP  
C  
      IF(ISTRAN(4).GE.1)THEN  
        PRINT *,'READING SFBSER.INP'
        OPEN(1,FILE='SFBSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,18  
          READ(1,*)  
        ENDDO  
        READ(1,*,IOSTAT=ISO) MSFSER,TCSFSER,TASFSER,TSRSF,TSSSF,  
     &      ISSFLDN,ISSFLFE,SFLKILL  
        IF(ISO.GT.0) GOTO 960  
        DO M=1,MSFSER  
          READ(1,*,IOSTAT=ISO) TSFSER(M),RKDSFL(M),WSFLST(M),WSFLSM(M),  
     &        DSFLMN(M),DSFLMX(M),SFNTBE(M),SFATBT(M)  
          IF(ISO.GT.0) GOTO 960  
        ENDDO  
        CLOSE(1)  
      ENDIF  
C  
C **  READ VEGETATION DATA FROM VEGE.INP AND VEGSER.INP  
C  
      IF(ISVEG.GE.1)THEN
        PRINT *,'READING VEGE.INP'
        OPEN(1,FILE='VEGE.INP',STATUS='UNKNOWN')  
        DO NS=1,12  
          READ(1,*)  
        ENDDO  
        READ(1,*)MVEGTYP,MVEGOW,NVEGSER,UVEGSCL
        MCOUNT=0 !Reset macroalgae count
        DO M=1,MVEGTYP  
         READ(1,*,ERR=3120)IDUM,NVEGSERV(M),RDLPSQ(M),VEGDIAMTMP,VEGHPTMP,  
     &        ALPVEG(M),BETVEG(M),GAMVEG(M),SCVEG(M)
!SCJ If macroalgae are to be in more than just the bottom layer, set GAMVEG equal to the elevation where growth starts (bottom of kelp farm).
!SCJ If macroalgae should actually start at 0 elevation, then set GAMVEG = -1.0E-6.
          IF(GAMVEG(M)==0.0)THEN !Flag for regular vegetation
            BPVEG(M)=VEGDIAMTMP
            HPVEG(M)=VEGHPTMP
            BDLTMP=BPVEG(M)*BPVEG(M)*RDLPSQ(M)  
            PVEGX(M)=1.-BETVEG(M)*BDLTMP  !Never used  
            PVEGY(M)=1.-BETVEG(M)*BDLTMP  !Never used
            PVEGZ(M)=1.-ALPVEG(M)*BDLTMP
            BDLPSQ(M)=BPVEG(M)*RDLPSQ(M)
          ELSE !Macroalgae
            DO L=2,LA
              IF(M/=MVEGL(L))CYCLE
              MACDIAM(L)=VEGDIAMTMP !Kelp diameter changes with growth
              ZMAXMAC(L)=VEGHPTMP !Kelp length changes with growth
              ZMINMAC(L)=GAMVEG(M) !Bottom of kelp farm doesn't change, but L dimension kept for consistency (unused GAMVEG variable used for this)
              IF(ZMAXMAC(L)<=ZMINMAC(L))THEN
                  PRINT*,'MacMax (HPVEG)<MacMin (GAMVEG) in VEGE.INP',L
                  PAUSE
              ENDIF
              BDLTMP=MACDIAM(L)*MACDIAM(L)*RDLPSQ(M)!Normalized area of MACROALGAE in this cell (-)
              IF(ALPVEG(M)*BDLTMP>1.0)THEN
                  PRINT*,'Macroalgae (including blocking factor) area exceeds cell area. Check VEGE.INP: RDLPSQ, BPVEG, ALPVEG ',PAUSE
              ENDIF
              MVEGZ(L)=MAX(0.0,1.0-ALPVEG(M)*BDLTMP)!Blockage ratio in Z (-)
              MACAD(L)=RDLPSQ(M)*MACDIAM(L)        !a*d (m) used to calculate drag, see papers written on this James and O'Donncha EFM 2019 BDLPSQ
              MCOUNT=MCOUNT+1
              IJLMAC(MCOUNT,1)=IL(L) !Not used
              IJLMAC(MCOUNT,2)=JL(L) !Not used
              IJLMAC(MCOUNT,3)=L
            ENDDO
          ENDIF
        ENDDO  
        CLOSE(1)  
        IF(NVEGSER.GT.0)THEN  
          DO M=1,NVEGSER  
            MVEGTLAST(M)=1  
          ENDDO  
        ENDIF
      ENDIF
	!!!BEGIN SCJ BLOCK
      IF(LMHK)THEN !MHK devices
        PRINT *,'READING MHK.INP'
        OPEN(1,FILE='MHK.INP',STATUS='UNKNOWN')  
        DO NS=1,29  
          READ(1,*)  
        ENDDO  
        READ(1,*,ERR=3122)MHKTYP,NFLAGPWR,UPSTREAM,OUTPUTFLAG
        IF(NFLAGPWR==1)THEN
          DO M=1,MHKTYP  
            READ(1,*,ERR=3122)WIDTHMHK(M),WIDTHSUP(M),
     &BOFFMHK(M),BOFFSUP(M),TOFFMHK(M),TOFFSUP(M),
     &CTMHK(M),CDSUP(M),VMINCUT(M),VMAXCUT(M),DENMHK(M)
            CTMHK(M)=CTMHK(M)*DENMHK(M)
            CDSUP(M)=CDSUP(M)*DENMHK(M)
		  DO L=2,LA
		    IF(M+90==MVEGL(L))THEN      ! Can have number of different MHK types (dimensions etc.) here we have 5 MHK devices but all are equivalent, so only need to scan across cells containing MHK devices with the same MHK parameters
		      ZMINMHK(M,L)=BELV(L)+BOFFMHK(M)
		      ZMAXMHK(M,L)=BELV(L)+TOFFMHK(M)
		      ZMINSUP(M,L)=BELV(L)+BOFFSUP(M)
		      ZMAXSUP(M,L)=BELV(L)+TOFFSUP(M)
	          DIAMMHK=ZMAXMHK(M,L)-ZMINMHK(M,L)
                IF(DIAMMHK<0.0)THEN
		        PRINT*,'MHK ZMIN > ZMAX'
	            STOP
	          ENDIF
	        ENDIF
		  ENDDO
          ENDDO
          READ(1,*) !skip the header line
          READ(1,*)BETAMHK_D,BETAMHK_P,CE4MHK,PB_COEF
        ELSEIF(NFLAGPWR==2)THEN
          PRINT*,'Not available yet'
        ELSEIF(NFLAGPWR==3)THEN
          DO M=1,MHKTYP
            READ(1,*,ERR=3122)WIDTHMHK(M),WIDTHSUP(M),
     &BOFFMHK(M),HEIGHTMHK(M),HEIGHTSUP(M),REFELEV(M),
     &CTMHK(M),CDSUP(M),VMINCUT(M),VMAXCUT(M),DENMHK(M)
		  CTMHK(M)=CTMHK(M)*DENMHK(M)
		  CDSUP(M)=CDSUP(M)*DENMHK(M)
		  DO L=2,LA
		    IF(M+90==MVEGL(L))THEN 
		      ZMINMHK(M,L)=BELV(L)+REFELEV(M)
		      ZMAXMHK(M,L)=BELV(L)+REFELEV(M)+HEIGHTMHK(M)
		      ZMINSUP(M,L)=BELV(L)
		      ZMAXSUP(M,L)=BELV(L)+REFELEV(M)+HEIGHTSUP(M)
	          DIAMMHK=ZMAXMHK(M,L)-ZMINMHK(M,L)
                IF(DIAMMHK<0.0)THEN
		        PRINT*,'MHK ZMIN > ZMAX'
	            STOP
	          ENDIF
	        ENDIF
		  ENDDO
            READ(1,*) !skip the header line
            READ(1,*)BETAMHK_D,BETAMHK_P,CE4MHK,PB_COEF
          ENDDO
        ENDIF
        CLOSE(1)    
        DO L=2,LA
	    IF(MVEGL(L)>90)THEN
	      IF((DIAMMHK>DXP(L).OR.DIAMMHK>DYP(L))
     &.AND.DENMHK(MVEGL(L)-90)>1.0)THEN
	        PRINT*,'MHK DIAMETER EXCEEDS CELL SIZE'
	        PRINT*,'AND DENSITY >= 1'
	        STOP
	      ENDIF
	    ENDIF
        ENDDO 
      ENDIF
	!!!END SCJ BLOCK
      GOTO 3124  
 3120 WRITE(6,3121)  
 3121 FORMAT('  READ ERROR FOR FILE VEGE.INP ')  
      STOP
 3122 WRITE(6,3123)  
 3123 FORMAT('  READ ERROR FOR FILE MHK.INP ')  
      STOP 
 3124 CONTINUE  
C_____________________________________________________________________
C BEGIN FEARGHAL O DONNCHA AQUACULTURE COMPONENT
C READ DATA RELATED TO aquaculture installations
C______________________________________________________________________
! check if input aquaculture configuration file exist
      INQUIRE(FILE='INPUTAQUA.INP',exist = exist)
      IF (exist) then ! open the aquaculture input file
         OPEN(100,FILE='INPUTAQUA.INP',STATUS='OLD')
         do i = 1,9
            read(100,*)   ! skip header lines
         end do
         read (100,*) Aquaincl,aquadiam,aquaden,aquadrag,haqua
         do i = 1,7
            read(100,*) ! skip header lines
         end do
         read(100,*) Aquaturb,BETADRP_P,BETADRP_D,CE4DRP ! turbulence characteristics
         close(100)
! aquaincl = flag for the inclusion of aquacultur
! thetapar = area reduction coefficient (=1-DropperArea*DropperDensity)
! aquadrag = aquaculture instl. cylinder drag coefficent (0.8 - 2.0)
! aquadiam = aquaculture instl. diam.
! aquaden = number of aquaculture canopy droppers per sq. mtr.
! kcdrp = modification applied from layer KCDRP  to upper layer
! aquaturb = flag to implement aquaculture turbulence effects
! ctkdrop = efficiency of the production of teh energy of turbulence
! ctldrop = effficiency of dissipation of turbulent length scale
! haqua = height of aquaculture installation (distance from surface
         CLOSE(100)
      END IF       

      IF (aquaincl == 1 ) THEN 
         call structident
         WRITE(*,*) "AQUACULTURE INITIALIZATION DONE"
      END IF
!*********************************************************************
! END Aquaculture related call

      IF(NVEGSER.GE.1)THEN  
        PRINT *,'READING VEGSER.INP'
        OPEN(1,FILE='VEGSER.INP',STATUS='UNKNOWN')  
C  
C **  SKIP OVER TITLE AND AND HEADER LINES  
C  
        DO IS=1,8  
          READ(1,*)  
        ENDDO  
        DO NS=1,NVEGSER  
          READ(1,*,IOSTAT=ISO) MVEGSER(NS),TCVEGSER(NS),TAVEGSER(NS)  
          IF(ISO.GT.0) GOTO 7120  
          DO M=1,MVEGSER(NS)  
            READ(1,*,IOSTAT=ISO)TVEGSER(M,NS),VEGSERR(M,NS),  
     &          VEGSERB(M,NS),VEGSERH(M,NS)  
            IF(ISO.GT.0) GOTO 7120  
            TVEGSER(M,NS)=TVEGSER(M,NS)+TAVEGSER(NS)  
          ENDDO  
        ENDDO  
        CLOSE(1)  
C  
C **  REINITIALIZE CLASSES HAVING TIME SERIES INFORMATION  
C  
        DO M=1,MVEGTYP  
          IF(NVEGSERV(M).GT.0)THEN  
            NS=NVEGSERV(M)  
            RDLPSQ(M)=VEGSERR(1,NS)  
            BPVEG(M)=VEGSERB(1,NS)  
            HPVEG(M)=VEGSERH(1,NS)  
            BDLTMP=BPVEG(M)*BPVEG(M)*RDLPSQ(M)  
            PVEGX(M)=1.-BETVEG(M)*BDLTMP  
            PVEGY(M)=1.-BETVEG(M)*BDLTMP  
            PVEGZ(M)=1.-ALPVEG(M)*BDLTMP  
            BDLPSQ(M)=BPVEG(M)*RDLPSQ(M)  
          ENDIF  
        ENDDO  
      ENDIF  
      INQUIRE(FILE='PERTURBATIONS.INP',exist = exist)
      IF (exist) then
        OPEN(1,FILE='PERTURBATIONS.INP',STATUS='OLD')
        DO I = 1,10
          READ(1,*)
        END DO
        READ(1,*) BASE,PERT_MIN,PERT_MAX
      END IF
      CLOSE(1)
  
      IF (ISFAVB == 4) THEN   ! READ TURBULENCE CALIBRATION PARAMETERS
                              ! FOR VINCON-LEITE IMPLEMENTATION
         OPEN(1,FILE='VINCLETPAR.INP',STATUS='OLD')
         DO I = 1,13
           READ(1,*) 
         END DO
         READ(1,*,IOSTAT=ISO) CP1, CP2, CP3, ALPHA_VLET, P2_MET, P2_HYPO
         IF(ISO.GT.0) GOTO 980
         CLOSE(1)
      END IF
! open input tidal turbine file

C
C ** CHECK IF DATA ASSIMILATION IS IMPLEMENTED AND READ INPUT DATA
C
! compiler flags specified in makefile overwrite input file flags (cannot invoke without LA)
#ifdef key_da
      INQUIRE(FILE='DA.INP',exist = exist)
      IF (exist) THEN
        OPEN(100,FILE='DA.INP',STATUS='OLD')
        DO I =1,11
          READ(100,*) ! SKIP HEADER FILES
        END DO
        READ(100,*) STRING_DESCR, IDA_FLAG
        write(*,*) 'IDAFLAG =', IDA_FLAG
        IF (IDA_FLAG ==1) THEN
          write(*,*) 'DATA ASSIMILATION INVOKED -- READ CONFIGURATIONS FROM FILE DA.INP'
          READ(100, *)  ! Skip 2 lines
          READ(100, *)
          READ(100,*) STRING_DESCR, IBEG_DA    !   \
          READ(100,*) STRING_DESCR, IEND_DA    !    \ SPATIAL EXTENTS OF
          READ(100,*) STRING_DESCR, JBEG_DA    !    / THE DA ROUTINE
          READ(100,*) STRING_DESCR, JEND_DA    !   /  (RECTANGULAR GRID)
          READ(100,*) STRING_DESCR, NDAPOINTS  !  NUMBER OF EFDC GRID POINTS TO ASSIMILATE
          READ(100,*) STRING_DESCR, DAFREQ

          write(*,*) 'DAFREQ = ', DAFREQ
C **  BASED ON DAFREQ COMPUTE TIMESTEP FREQUENCY TO CALL DATA ASSIMILATION
          NDATA_ASSIM = DAFREQ * 3600/DT  ! DA FREQUENCY PER HOUR
          READ(100,*)
          READ(100,*) STRING_DESCR, PMATRIX_R1
          READ(100,*) STRING_DESCR, PMATRIX_R2
          READ(100,*) STRING_DESCR, PMATRIX_A
          READ(100,*) STRING_DESCR, EKPROJ
          ALLOCATE(IBLUE(NDAPOINTS))
          ALLOCATE(JBLUE(NDAPOINTS))
          ALLOCATE(LBLUE(NDAPOINTS))
        END IF
      END IF
#endif
!
 
      GOTO 7122  
 7120 WRITE(6,7121)  
 7121 FORMAT('  READ ERROR FOR FILE VEGSER.INP ')  
      STOP  
 7122 CONTINUE  
      GOTO 3000  
  
C **  WRITE READ ERROR FOR OTHER INPUT FILES AND TERMINATE RUN  
C  
  800 WRITE(6,801)  
      WRITE(8,801)  
  801 FORMAT('  READ ERROR FOR FILE CELL.INP ')  
      STOP  
  820 WRITE(6,821)  
      WRITE(8,821)  
  821 FORMAT('  READ ERROR FOR FILE DEPTH.INP ')  
      STOP  
  830 WRITE(6,831)  
      WRITE(8,831)  
  831 FORMAT('  READ ERROR FOR FILE DXDY.INP ')  
      STOP  
  840 WRITE(6,841)  
      WRITE(8,841)  
  841 FORMAT('  READ ERROR FOR FILE SALT.INP ')  
      STOP  
  842 WRITE(6,843)  
      WRITE(8,843)  
  843 FORMAT('  READ ERROR FOR FILE TEMP.INP ')  
      STOP  
  844 WRITE(6,845)  
      WRITE(8,845)  
  845 FORMAT('  READ ERROR FOR FILE DYE.INP ')  
      STOP  
  846 WRITE(6,847)  
      WRITE(8,847)  
  847 FORMAT('  READ ERROR FOR FILE SFL.INP ')  
      STOP  
  848 WRITE(6,849)  
      WRITE(8,849)  
  849 FORMAT('  READ ERROR FOR FILE TOXW.INP ')  
      STOP  
  852 WRITE(6,853)  
      WRITE(8,853)  
  853 FORMAT('  READ ERROR FOR FILE TOXB.INP ')  
      STOP  
  850 WRITE(6,851)  
      WRITE(8,851)  
  851 FORMAT('  READ ERROR FOR FILE PSER.INP OR UVSER.INP ')  
      STOP  
  854 WRITE(6,855)  
      WRITE(8,855)  
  855 FORMAT('  READ ERROR FOR FILE SEDW.INP ')  
      STOP  
  856 WRITE(6,857)  
      WRITE(8,857)  
  857 FORMAT('  READ ERROR FOR FILE SEDB.INP ')  
      STOP  
  858 WRITE(6,859)  
      WRITE(8,859)  
  859 FORMAT('  READ ERROR FOR FILE SNDW.INP ')  
      STOP  
  862 WRITE(6,863)  
      WRITE(8,863)  
  863 FORMAT('  READ ERROR FOR FILE SNDB.INP ')  
      STOP  
  860 WRITE(6,861)  
      WRITE(8,861)  
  861 FORMAT('  READ ERROR FOR FILE QSER.INP ')  
      STOP  
  865 WRITE(6,866)  
      WRITE(8,866)  
  866 FORMAT('  READ ERROR FOR FILE QWRS.INP ')  
      STOP  
  870 WRITE(6,871)  
      WRITE(8,871)  
  871 FORMAT('  READ ERROR FOR FILE SSER.INP ')  
      STOP  
  880 WRITE(6,881)  
      WRITE(8,881)  
  881 FORMAT('  READ ERROR FOR FILE TSER.INP ')  
      STOP  
  890 WRITE(6,891)  
      WRITE(8,891)  
  891 FORMAT('  READ ERROR FOR FILE DSER.INP ')  
      STOP  
  900 WRITE(6,901)  
      WRITE(8,901)  
  901 FORMAT('  READ ERROR FOR FILE SDSER.INP ')  
      STOP  
  902 WRITE(6,903)  
      WRITE(8,903)  
  903 FORMAT('  READ ERROR FOR FILE SNSER.INP ')  
      STOP  
  904 WRITE(6,905)  
      WRITE(8,905)  
  905 FORMAT('  READ ERROR FOR FILE TXSER.INP ')  
      STOP  
  910 WRITE(6,911)  
      WRITE(8,911)  
  911 FORMAT('  READ ERROR FOR FILE SFSER.INP ')  
      STOP  
  920 WRITE(6,921)  
      WRITE(8,921)  
  921 FORMAT('  READ ERROR FOR FILE QCTL.INP ')  
      STOP  
  940 WRITE(6,941)  
      WRITE(8,941)  
  941 FORMAT('  READ ERROR FOR FILE ASER.INP ')  
      STOP  
  950 WRITE(6,951)  
      WRITE(8,951)  
  951 FORMAT('  READ ERROR FOR FILE MAPPGNS.INP ')  
      STOP  
  960 WRITE(6,961)  
      WRITE(8,961)  
  961 FORMAT('  READ ERROR FOR FILE SFBSER.INP ')  
      STOP  
  970 WRITE(6,971)  
      WRITE(8,971)  
  971 FORMAT('  READ ERROR FOR FILE TIDASM.INP ')  
      STOP  
  980 WRITE(6,981)  
      WRITE(8,981)  
  981 FORMAT('  READ ERROR FOR FILE VINCLETPAR.INP ')  
      STOP  
 3000 CONTINUE  
      RETURN  
      END  

      ! *** DSLLC UTIL
      FUNCTION PARSE_REAL(INLINE)

        CHARACTER*(*) INLINE
        CHARACTER*15  CVAL
        
        ILEN=LEN_TRIM(INLINE)
        PARSE_REAL=0.
        DO IC=1,ILEN
          IF(INLINE(IC:IC).EQ.':')THEN
            DO IPOS=IC+1,ILEN
              IF(INLINE(IPOS:IPOS).NE.' ')EXIT
            ENDDO
            IF(IPOS.GT.ILEN)RETURN
            CVAL=INLINE(IPOS:ILEN)
            ILEN2=LEN_TRIM(CVAL)
            DO JC=1,ILEN2
              IF(CVAL(JC:JC).EQ.' '.OR.CVAL(JC:JC).EQ.','
     &             .OR.JC.EQ.ILEN2)THEN
                READ(CVAL(1:JC),'(F12.1)',ERR=999)PARSE_REAL
                RETURN
              ENDIF
            ENDDO
          ENDIF
        ENDDO
  
  999   print *, ' error parsing real'
        RETURN
      END

      FUNCTION PARSE_LOGICAL(INLINE)

        CHARACTER*(*) INLINE
        CHARACTER*12  CVAL
        LOGICAL       PARSE_LOGICAL
        
        ILEN=LEN_TRIM(INLINE)
        DO IC=1,ILEN
          IF(INLINE(IC:IC).EQ.':')THEN
            DO IPOS=IC+1,ILEN
              IF(INLINE(IPOS:IPOS).NE.' ')EXIT
            ENDDO
            IF(IPOS.GT.ILEN)RETURN
            CVAL=INLINE(IPOS:ILEN)
            ILEN2=LEN_TRIM(CVAL)
            DO JC=1,ILEN2
              IF(CVAL(JC:JC).EQ.' '.OR.CVAL(JC:JC).EQ.','
     &             .OR.JC.EQ.ILEN2)THEN
                IF(CVAL(1:1).EQ.'T'.OR.CVAL(1:1).EQ.'Y')THEN
                  PARSE_LOGICAL=.TRUE.
                  RETURN
                ENDIF
              ENDIF
            ENDDO
          ENDIF
        ENDDO
        PARSE_LOGICAL=.FALSE.
  900   FORMAT(L1)  
  999   RETURN
      END

      SUBROUTINE CHECK_OUTPUT_FREQUENCY
      USE GLOBAL
      IMPLICIT NONE

      ! To allow effective mapping of variables to the Netcdf file
      ! we need consistency between output frequency for
      ! velocity, surface elevation and transported constitutents.
      ! These are all written to a single file efdcout_YYYYMMDDHHMMSS.nc
      ! at equivalent time frequency. Hence need to synchronise
      ! NPPPH (elevation), NPVPH(velocity) and NPSPH(scalar transport)
      ! Set write frequency to maximum of user specified values and print
      ! warning to screen
      INTEGER NPSPH_MAX

      ! 1) Find maximum write frequency of scalar transport
      NPSPH_MAX = 0
      DO N=1,7
        IF (ISSPH(N) > 0) NPSPH_MAX = MAX(NPSPH_MAX, NPSPH(N))
      END DO
      ! We need to allow for situation where user outputs solute transport or
      ! surface elevation but not velocities
      IF (ISVPH >0 .AND. ISPPH > 0) THEN
        IF (NPVPH /= NPPPH) THEN
          IF (PARTID == MASTER_TASK ) THEN
            PRINT *, ' USER DEFINED WRITE INTERVAL FOR ELEVATION AND VELOCITY NOT EQUAL '
            PRINT *, ' SPECIFYING WRITE FREQUENCY AS MAXIMUM VALUE = ', MAX(NPVPH,NPPPH)
          END IF
          NPVPH = MAX(NPVPH,NPPPH)
        END IF
      END IF
      IF (ISPPH > 0) NPVPH = MAX(NPVPH,NPPPH)
      IF (NPSPH_MAX > 0) THEN ! (i.e. write if solute transport is active)
        IF (NPVPH /= NPSPH_MAX) THEN
          IF (PARTID == MASTER_TASK ) THEN
            PRINT *, ' USER DEFINED WRITE INTERVAL FORTRANSPORTED VARIABLES AND VEL/EL NOT EQUAL '
            PRINT *, ' SPECIFYING WRITE FREQUENCY AS MAXIMUM VALUE = ',MAX(NPVPH,NPSPH_MAX)
          END IF
          NPVPH = MAX(NPVPH,NPSPH_MAX)
        END IF
      END IF
      IF (PARTID == MASTER_TASK ) THEN
        PRINT *, ' HORIZONTAL MAPS OF USER DEFINED VARIABLES (VEL, ELEV, TRANSPORT)'
        PRINT *, ' WILL BE WRITTEN', NPVPH, '  TIMES PER REFERENCE PERIOD'
      END IF
      RETURN
      END SUBROUTINE CHECK_OUTPUT_FREQUENCY

